ProteoWizard
SpectrumList_MZRefinerTest.cpp
Go to the documentation of this file.
1 //
2 // $Id: SpectrumList_MZRefinerTest.cpp 8985 2015-10-12 22:38:57Z chambm $
3 //
4 //
5 // Original author: Bryson Gibbons <bryson.gibbons@pnnl.gov>
6 //
7 // Copyright 2014 Pacific Northwest National Laboratory
8 // Richland, WA 99352
9 //
10 // Licensed under the Apache License, Version 2.0 (the "License");
11 // you may not use this file except in compliance with the License.
12 // You may obtain a copy of the License at
13 //
14 // http://www.apache.org/licenses/LICENSE-2.0
15 //
16 // Unless required by applicable law or agreed to in writing, software
17 // distributed under the License is distributed on an "AS IS" BASIS,
18 // WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
19 // See the License for the specific language governing permissions and
20 // limitations under the License.
21 //
22 
23 
28 #include "boost/filesystem/path.hpp"
30 #include <cstring>
31 
32 
33 using namespace pwiz::cv;
34 using namespace pwiz::msdata;
35 using namespace pwiz::util;
36 using namespace pwiz::analysis;
37 namespace bfs = boost::filesystem;
38 
39 
40 ostream* os_ = 0;
41 
42 // Check scan metadata and a small sample of the m/z data for high-res scans.
43 void verifyScanInfo(const Spectrum& spectrum, const double& epsilon, double basePeakMZ, double lowestObservedMZ, double highestObservedMZ, int mzArrayIndex1, double mzArrayValue1, int mzArrayIndex2, double mzArrayValue2)
44 {
45  unit_assert(spectrum.hasBinaryData());
46  const BinaryDataArrayPtr binaryData = spectrum.getMZArray();
47 
48  if (os_)
49  {
50  *os_ << "[verifyScanInfo] " << spectrum.index << " " << spectrum.id << " "
51  << basePeakMZ << " " << lowestObservedMZ << " " << highestObservedMZ << " "
52  << mzArrayValue1 << " " << mzArrayValue2 << ": "
53  << spectrum.cvParam(MS_base_peak_m_z).value << " " << spectrum.cvParam(MS_lowest_observed_m_z).value << " "
54  << spectrum.cvParam(MS_highest_observed_m_z).value << " " << binaryData->data[mzArrayIndex1] << " "
55  << binaryData->data[mzArrayIndex2] << endl;
56  }
57 
58  unit_assert_equal(spectrum.cvParam(MS_base_peak_m_z).valueAs<double>(), basePeakMZ, epsilon);
59  unit_assert_equal(spectrum.cvParam(MS_lowest_observed_m_z).valueAs<double>(), lowestObservedMZ, epsilon);
60  unit_assert_equal(spectrum.cvParam(MS_highest_observed_m_z).valueAs<double>(), highestObservedMZ, epsilon);
61  unit_assert_equal(binaryData->data[mzArrayIndex1], mzArrayValue1, epsilon);
62  unit_assert_equal(binaryData->data[mzArrayIndex2], mzArrayValue2, epsilon);
63 }
64 
65 // Check scan precursor metadata for MS/MS scans
66 void verifyPrecursorInfo(const Spectrum& spectrum, const double& epsilon, double precursorMZ, double isolationWindowTarget)
67 {
68  unit_assert(!spectrum.precursors.empty());
69  const Precursor& precursor = spectrum.precursors[0];
70  unit_assert(!precursor.selectedIons.empty());
71  const SelectedIon& selectedIon = precursor.selectedIons[0];
72  unit_assert(!precursor.isolationWindow.empty());
73  const IsolationWindow& isoWindow = precursor.isolationWindow;
74 
75  if (os_)
76  {
77  *os_ << "[verifyPrecursorInfo] " << spectrum.index << " " << spectrum.id << " "
78  << precursorMZ << " " << isolationWindowTarget << ": "
79  << selectedIon.cvParam(MS_selected_ion_m_z).value << " " << isoWindow.cvParam(MS_isolation_window_target_m_z).value << endl;
80  }
81 
82  unit_assert_equal(selectedIon.cvParam(MS_selected_ion_m_z).valueAs<double>(), precursorMZ, epsilon);
83  unit_assert_equal(isoWindow.cvParam(MS_isolation_window_target_m_z).valueAs<double>(), isolationWindowTarget, epsilon);
84 }
85 
86 void testShift(const bfs::path& datadir)
87 {
88  MSDataFile msd((datadir / "JD_06232014_sample4_C.mzML").string());
89 
90  unit_assert(msd.run.spectrumListPtr.get() && msd.run.spectrumListPtr->size()==610);
91  if (os_) *os_ << "original spectra:\n";
92  // Provided mzML file is high-res/high-res
93  double epsilon = 1e-4;
94  // MS1 scans 0, 224, 398 (0 and 224
95  verifyScanInfo(*msd.run.spectrumListPtr->spectrum(0, true), epsilon, 371.09958, 300.14306, 1568.55126, 30, 303.64633, 1200, 416.24838);
96  verifyScanInfo(*msd.run.spectrumListPtr->spectrum(224, true), epsilon, 558.30688, 301.05908, 1522.72473, 200, 407.26425, 1500, 724.32824);
97  verifyScanInfo(*msd.run.spectrumListPtr->spectrum(10, true), epsilon, 530.32782, 74.06039, 887.42852, 41, 188.11117, 93, 442.22839);
98  verifyPrecursorInfo(*msd.run.spectrumListPtr->spectrum(10), epsilon, 530.26684, 530.27);
99  verifyScanInfo(*msd.run.spectrumListPtr->spectrum(173, true), epsilon, 141.10162, 87.05542, 1187.53137, 63, 248.15817, 116, 887.44793);
100  verifyPrecursorInfo(*msd.run.spectrumListPtr->spectrum(173), epsilon, 629.30160, 629.3);
101  verifyScanInfo(*msd.run.spectrumListPtr->spectrum(346, true), epsilon, 848.45895, 116.00368, 1454.73327, 16, 185.16418, 95, 862.43109);
102  verifyPrecursorInfo(*msd.run.spectrumListPtr->spectrum(346), epsilon, 840.45480, 840.45);
103  verifyScanInfo(*msd.run.spectrumListPtr->spectrum(470, true), epsilon, 249.15857, 119.04895, 1402.77331, 23, 217.08113, 102, 1154.59863);
104  verifyPrecursorInfo(*msd.run.spectrumListPtr->spectrum(470), epsilon, 838.96706, 838.97);
105  verifyScanInfo(*msd.run.spectrumListPtr->spectrum(551, true), epsilon, 1048.55047, 155.08105, 1321.67761, 50, 368.19134, 104, 941.96954);
106  verifyPrecursorInfo(*msd.run.spectrumListPtr->spectrum(551), epsilon, 739.69935, 740.03);
107 
108  shared_ptr<SpectrumList_MZRefiner> spectrumListMZRefined(
109  new SpectrumList_MZRefiner(msd, (datadir / "JD_06232014_sample4_C.mzid").string(), "specEValue", "-1e-10", IntegerSet(1, 2)));
110 
111  unit_assert(spectrumListMZRefined->size() == 610);
112  if (os_) *os_ << "refined spectra:\n";
113  epsilon = 1e-2; // Increase the tolerance a little bit for the refined results (add some forgiveness with algorithm updates...)
114  verifyScanInfo(*msd.run.spectrumListPtr->spectrum(0, true), epsilon, 371.10060, 300.14388, 1568.55631, 30, 303.64715, 1200, 416.24951);
115  verifyScanInfo(*msd.run.spectrumListPtr->spectrum(224, true), epsilon, 558.30841, 301.05990, 1522.72962, 200, 407.26538, 1500, 724.33007);
116  verifyScanInfo(*msd.run.spectrumListPtr->spectrum(10, true), epsilon, 530.32928, 74.06059, 887.43126, 41, 188.11169, 93, 442.22961);
117  verifyPrecursorInfo(*msd.run.spectrumListPtr->spectrum(10), epsilon, 530.26830, 530.27145);
118  verifyScanInfo(*msd.run.spectrumListPtr->spectrum(173, true), epsilon, 141.10200, 87.05566, 1187.53519, 63, 248.15885, 116, 887.45068);
119  verifyPrecursorInfo(*msd.run.spectrumListPtr->spectrum(173), epsilon, 629.30333, 629.30172);
120  verifyScanInfo(*msd.run.spectrumListPtr->spectrum(346, true), epsilon, 848.46155, 116.00400, 1454.73795, 16, 185.16468, 95, 862.43371);
121  verifyPrecursorInfo(*msd.run.spectrumListPtr->spectrum(346), epsilon, 840.45738, 840.45257);
122  verifyScanInfo(*msd.run.spectrumListPtr->spectrum(470, true), epsilon, 249.15926, 119.04927, 1402.77782, 23, 217.08172, 102, 1154.60229);
123  verifyPrecursorInfo(*msd.run.spectrumListPtr->spectrum(470), epsilon, 838.96963, 838.97257);
124  verifyScanInfo(*msd.run.spectrumListPtr->spectrum(551, true), epsilon, 1048.55384, 155.08147, 1321.68186, 50, 368.19235, 104, 941.97253);
125  verifyPrecursorInfo(*msd.run.spectrumListPtr->spectrum(551), epsilon, 739.70141, 740.03206);
126 }
127 
128 
129 void test(const bfs::path& datadir)
130 {
131  testShift(datadir);
132 }
133 
134 
135 int main(int argc, char* argv[])
136 {
137  TEST_PROLOG(argc, argv)
138 
139  try
140  {
141  bfs::path datadir = ".";
142 
143  // grab the parent directory for the test files.
144  for (int i=1; i<argc; i++)
145  {
146  if (!strcmp(argv[i],"-v"))
147  os_ = &cout;
148  else
149  // hack to allow running unit test from a different directory:
150  // Jamfile passes full path to specified input file.
151  // we want the path, so we can ignore filename
152  datadir = bfs::path(argv[i]).branch_path();
153  }
154 
155  if (argc>1 && !strcmp(argv[1],"-v")) os_ = &cout;
156  test(datadir);
157  }
158  catch (exception& e)
159  {
160  TEST_FAILED(e.what())
161  }
162  catch (...)
163  {
164  TEST_FAILED("Caught unknown exception.")
165  }
166 
168 }
169 
170 
std::string value
Definition: ParamTypes.hpp:47
CVParam cvParam(CVID cvid) const
finds cvid in the container:
a virtual container of integers, accessible via an iterator interface, stored as union of intervals ...
Definition: IntegerSet.hpp:37
The method of precursor ion selection and activation.
Definition: MSData.hpp:310
MS_highest_observed_m_z
highest observed m/z: Highest m/z value observed in the m/z array.
Definition: cv.hpp:2035
const double epsilon
Definition: DiffTest.cpp:41
ostream * os_
void verifyScanInfo(const Spectrum &spectrum, const double &epsilon, double basePeakMZ, double lowestObservedMZ, double highestObservedMZ, int mzArrayIndex1, double mzArrayValue1, int mzArrayIndex2, double mzArrayValue2)
std::vector< Precursor > precursors
list and descriptions of precursors to the spectrum currently being described.
Definition: MSData.hpp:519
#define TEST_EPILOG
Definition: unit.hpp:182
#define unit_assert_equal(x, y, epsilon)
Definition: unit.hpp:99
MS_selected_ion_m_z
selected ion m/z: Mass-to-charge ratio of an selected ion.
Definition: cv.hpp:2749
std::string id
a unique identifier for this spectrum. It should be expected that external files may use this identif...
Definition: MSData.hpp:475
MS_base_peak_m_z
base peak m/z: M/z value of the signal of highest intensity in the mass spectrum. ...
Definition: cv.hpp:1966
Run run
a run in mzML should correspond to a single, consecutive and coherent set of scans on an instrument...
Definition: MSData.hpp:882
bool hasBinaryData() const
returns true iff has nonnull and nonempty BinaryDataArrayPtr
Definition: MSData.hpp:534
size_t index
the zero-based, consecutive index of the spectrum in the SpectrumList.
Definition: MSData.hpp:472
MS_isolation_window_target_m_z
isolation window target m/z: The primary or reference m/z about which the isolation window is defined...
Definition: cv.hpp:3028
This element captures the isolation (or &#39;selection&#39;) window configured to isolate one or more precurs...
Definition: MSData.hpp:291
MSData object plus file I/O.
Definition: MSDataFile.hpp:40
void verifyPrecursorInfo(const Spectrum &spectrum, const double &epsilon, double precursorMZ, double isolationWindowTarget)
SpectrumListPtr spectrumListPtr
all mass spectra and the acquisitions underlying them are described and attached here. Subsidiary data arrays are also both described and attached here.
Definition: MSData.hpp:823
#define TEST_FAILED(x)
Definition: unit.hpp:176
BinaryDataArrayPtr getMZArray() const
get m/z array (may be null)
#define TEST_PROLOG(argc, argv)
Definition: unit.hpp:174
void testShift(const bfs::path &datadir)
The structure that captures the generation of a peak list (including the underlying acquisitions) ...
Definition: MSData.hpp:504
MS_lowest_observed_m_z
lowest observed m/z: Lowest m/z value observed in the m/z array.
Definition: cv.hpp:2038
value_type valueAs() const
templated value access with type conversion
Definition: ParamTypes.hpp:112
int main(int argc, char *argv[])
#define unit_assert(x)
Definition: unit.hpp:85
boost::shared_ptr< BinaryDataArray > BinaryDataArrayPtr
Definition: MSData.hpp:416
Definition: cv.hpp:91
void test(const bfs::path &datadir)