36 #ifndef OPENMS_FILTERING_CALIBRATION_INTERNALCALIBRATION_H 37 #define OPENMS_FILTERING_CALIBRATION_INTERNALCALIBRATION_H 84 template <
typename InputPeakType>
99 template <
typename InputPeakType>
115 template <
typename InputPeakType>
137 void calibrateMapGlobally(
const FeatureMap & feature_map,
FeatureMap & calibrated_feature_map, std::vector<PeptideIdentification> & ref_ids,
String trafo_file_name =
"");
141 template <
typename InputPeakType>
148 void makeLinearRegression_(std::vector<double> & observed_masses, std::vector<double> & theoretical_masses);
151 void checkReferenceIds_(std::vector<PeptideIdentification> & pep_ids);
154 void checkReferenceIds_(
const FeatureMap & feature_map);
164 template <
typename InputPeakType>
167 #ifdef DEBUG_CALIBRATION 172 std::cout <<
"Input is empty." << std::endl;
178 std::cout <<
"Attention: this function is assuming peak data." << std::endl;
180 calibrated_exp = exp;
182 Size num_ref_peaks = ref_masses.size();
183 bool use_ppm = (param_.getValue(
"mz_tolerance_unit") ==
"ppm") ?
true :
false;
184 double mz_tol = param_.getValue(
"mz_tolerance");
185 startProgress(0, exp.
size(),
"calibrate spectra");
187 for (
Size spec = 0; spec < exp.
size(); ++spec)
190 if (exp[spec].getMSLevel() != 1)
196 std::vector<double> corr_masses, rel_errors, found_ref_masses;
198 for (
Size peak = 0; peak < exp[spec].
size(); ++peak)
200 for (
Size ref_peak = 0; ref_peak < num_ref_peaks; ++ref_peak)
202 if (!use_ppm && fabs(exp[spec][peak].getMZ() - ref_masses[ref_peak]) < mz_tol)
204 found_ref_masses.push_back(ref_masses[ref_peak]);
205 corr_masses.push_back(exp[spec][peak].getMZ());
209 else if (use_ppm && fabs(exp[spec][peak].getMZ() - ref_masses[ref_peak]) / ref_masses[ref_peak] * 1e6 < mz_tol)
211 found_ref_masses.push_back(ref_masses[ref_peak]);
212 corr_masses.push_back(exp[spec][peak].getMZ());
220 std::cout <<
"spec: " << spec
221 <<
" less than 2 reference masses were detected within a reasonable error range\n";
222 std::cout <<
"This spectrum cannot be calibrated!\n";
227 for (
Size ref_peak = 0; ref_peak < found_ref_masses.size(); ++ref_peak)
229 rel_errors.push_back((found_ref_masses[ref_peak] - corr_masses[ref_peak]) / corr_masses[ref_peak] * 1e6);
232 makeLinearRegression_(corr_masses, found_ref_masses);
235 for (
unsigned int peak = 0; peak < calibrated_exp[spec].
size(); ++peak)
237 #ifdef DEBUG_CALIBRATION 238 std::cout << calibrated_exp[spec][peak].getMZ() <<
"\t";
240 double mz = calibrated_exp[spec][peak].getMZ();
241 mz = trafo_.apply(mz);
242 calibrated_exp[spec][peak].setMZ(mz);
243 #ifdef DEBUG_CALIBRATION 244 std::cout << calibrated_exp[spec][peak].getMZ() << std::endl;
253 template <
typename InputPeakType>
256 std::vector<PeptideIdentification> & ref_ids,
String trafo_file_name)
258 bool use_ppm = param_.getValue(
"mz_tolerance_unit") ==
"ppm" ?
true :
false;
259 double mz_tolerance = param_.getValue(
"mz_tolerance");
262 std::cout <<
"Input is empty." << std::endl;
268 std::cout <<
"Attention: this function is assuming peak data." << std::endl;
271 checkReferenceIds_(ref_ids);
273 std::vector<double> theoretical_masses, observed_masses;
274 for (
Size p_id = 0; p_id < ref_ids.size(); ++p_id)
276 for (
Size p_h = 0; p_h < ref_ids[p_id].getHits().size(); ++p_h)
278 Int charge = ref_ids[p_id].getHits()[p_h].getCharge();
279 double theo_mass = ref_ids[p_id].getHits()[p_h].getSequence().getMonoWeight(
Residue::Full, charge) / (
double)charge;
282 while (rt_iter != exp.
begin() && rt_iter->getMSLevel() != 1)
289 double dist = ref_ids[p_id].getMZ() - mz_iter->getMZ();
291 if ((mz_iter + 1) != rt_iter->
end()
292 && fabs((mz_iter + 1)->getMZ() - ref_ids[p_id].getMZ()) < fabs(dist)
293 && mz_iter != rt_iter->
begin()
294 && fabs((mz_iter - 1)->getMZ() - ref_ids[p_id].getMZ()) < fabs((mz_iter + 1)->getMZ() - ref_ids[p_id].getMZ()))
297 fabs((mz_iter + 1)->getMZ() - ref_ids[p_id].getMZ()) / ref_ids[p_id].getMZ() * 1e06 < mz_tolerance) ||
298 (!use_ppm && fabs((mz_iter + 1)->getMZ() - ref_ids[p_id].getMZ()) < mz_tolerance))
301 observed_masses.push_back((mz_iter + 1)->getMZ());
302 theoretical_masses.push_back(theo_mass);
307 else if (mz_iter != rt_iter->
begin()
308 && fabs((mz_iter - 1)->getMZ() - ref_ids[p_id].getMZ()) < fabs(dist))
311 fabs((mz_iter - 1)->getMZ() - ref_ids[p_id].getMZ()) / ref_ids[p_id].getMZ() * 1e06 < mz_tolerance) ||
312 (!use_ppm && fabs((mz_iter - 1)->getMZ() - ref_ids[p_id].getMZ()) < mz_tolerance))
315 observed_masses.push_back((mz_iter - 1)->getMZ());
316 theoretical_masses.push_back(theo_mass);
324 fabs((mz_iter)->getMZ() - ref_ids[p_id].getMZ()) / ref_ids[p_id].getMZ() * 1e06 < mz_tolerance) ||
325 (!use_ppm && fabs((mz_iter)->getMZ() - ref_ids[p_id].getMZ()) < mz_tolerance))
328 observed_masses.push_back(mz_iter->getMZ());
329 theoretical_masses.push_back(theo_mass);
337 makeLinearRegression_(observed_masses, theoretical_masses);
342 for (
Size spec = 0; spec < exp.
size(); ++spec)
345 if (exp[spec].getMSLevel() != 1)
347 calibrated_exp[spec] = exp[spec];
351 calibrated_exp[spec] = exp[spec];
353 for (
unsigned int peak = 0; peak < exp[spec].
size(); ++peak)
355 #ifdef DEBUG_CALIBRATION 356 std::cout << exp[spec][peak].getMZ() <<
"\t";
358 double mz = exp[spec][peak].getMZ();
359 mz = trafo_.apply(mz);
360 calibrated_exp[spec][peak].setMZ(mz);
361 #ifdef DEBUG_CALIBRATION 362 std::cout << calibrated_exp[spec][peak].getMZ() << std::endl;
367 if (trafo_file_name !=
"")
373 template <
typename InputPeakType>
378 std::cout <<
"Input is empty." << std::endl;
384 std::cout <<
"Attention: this function is assuming peak data." << std::endl;
388 Size num_ref_peaks = ref_masses.size();
389 bool use_ppm = (param_.getValue(
"mz_tolerance_unit") ==
"ppm") ?
true :
false;
390 double mz_tol = param_.getValue(
"mz_tolerance");
391 startProgress(0, exp.
size(),
"calibrate spectra");
392 std::vector<double> corr_masses, rel_errors, found_ref_masses;
395 for (
Size spec = 0; spec < exp.
size(); ++spec)
398 if (exp[spec].getMSLevel() != 1)
400 for (
Size peak = 0; peak < exp[spec].
size(); ++peak)
402 for (
Size ref_peak = 0; ref_peak < num_ref_peaks; ++ref_peak)
404 if (!use_ppm && fabs(exp[spec][peak].getMZ() - ref_masses[ref_peak]) < mz_tol)
406 found_ref_masses.push_back(ref_masses[ref_peak]);
407 corr_masses.push_back(exp[spec][peak].getMZ());
411 else if (use_ppm && fabs(exp[spec][peak].getMZ() - ref_masses[ref_peak]) / ref_masses[ref_peak] * 1e6 < mz_tol)
413 found_ref_masses.push_back(ref_masses[ref_peak]);
414 corr_masses.push_back(exp[spec][peak].getMZ());
423 std::cout <<
"Less than 2 reference masses were detected within a reasonable error range\n";
424 std::cout <<
"This spectrum cannot be calibrated!\n";
429 makeLinearRegression_(corr_masses, found_ref_masses);
434 for (
Size spec = 0; spec < exp.
size(); ++spec)
437 if (exp[spec].getMSLevel() != 1)
439 calibrated_exp[spec] = exp[spec];
444 calibrated_exp[spec] = exp[spec];
446 for (
unsigned int peak = 0; peak < exp[spec].
size(); ++peak)
448 #ifdef DEBUG_CALIBRATION 449 std::cout << exp[spec][peak].getMZ() <<
"\t";
451 double mz = exp[spec][peak].getMZ();
452 mz = trafo_.apply(mz);
453 calibrated_exp[spec][peak].setMZ(mz);
455 #ifdef DEBUG_CALIBRATION 456 std::cout << calibrated_exp[spec][peak].getMZ() << std::endl;
463 if (trafo_file_name !=
"")
471 #endif // OPENMS_FILTERING_CALIBRATION_INTERNALCALIBRATION_H
A more convenient string class.
Definition: String.h:57
Size size() const
Definition: MSExperiment.h:117
Peak data (also called centroided data or stick data)
Definition: SpectrumSettings.h:74
A container for features.
Definition: FeatureMap.h:93
A simple calibration method using linear interpolation of given reference masses. ...
Definition: InternalCalibration.h:61
ContainerType::const_iterator ConstIterator
Non-mutable iterator.
Definition: MSSpectrum.h:123
Iterator begin()
Definition: MSExperiment.h:147
void resize(Size s)
Definition: MSExperiment.h:122
Main OpenMS namespace.
Definition: FeatureDeconvolution.h:47
Iterator end()
Definition: MSExperiment.h:157
ConstIterator RTBegin(CoordinateType rt) const
Fast search for spectrum range begin.
Definition: MSExperiment.h:374
void calibrateMapGlobally(const MSExperiment< InputPeakType > &exp, MSExperiment< InputPeakType > &calibrated_exp, std::vector< double > &ref_masses, String trafo_file_name="")
Calibrate a peak map using given reference masses with one calibration function for the whole map...
Definition: InternalCalibration.h:374
~InternalCalibration()
Destructor.
Definition: InternalCalibration.h:70
Definition: Residue.h:361
In-Memory representation of a mass spectrometry experiment.
Definition: MSExperiment.h:69
void calibrateMapSpectrumwise(const MSExperiment< InputPeakType > &exp, MSExperiment< InputPeakType > &calibrated_exp, std::vector< double > &ref_masses)
Calibrate a peak map using given reference masses with a separate calibration function for each spect...
Definition: InternalCalibration.h:165
std::vector< SpectrumType >::const_iterator ConstIterator
Non-mutable iterator.
Definition: MSExperiment.h:103
Base class for all classes that want to report their progress.
Definition: ProgressLogger.h:55
bool empty() const
Definition: MSExperiment.h:127
Int writtenDigits< double >(const double &)
Number of digits commonly used for writing a double (a.k.a. precision).
Definition: Types.h:213
A base class for all classes handling default parameters.
Definition: DefaultParamHandler.h:92
int Int
Signed integer type.
Definition: Types.h:96
Description of the experimental settings.
Definition: ExperimentalSettings.h:59
TransformationDescription trafo_
here the transformation is stored
Definition: InternalCalibration.h:160