ProteoWizard
ChemistryTest.cpp
Go to the documentation of this file.
1 //
2 // $Id: ChemistryTest.cpp 6865 2014-10-31 21:47:12Z chambm $
3 //
4 //
5 // Original author: Darren Kessner <darren@proteowizard.org>
6 //
7 // Copyright 2006 Louis Warschaw Prostate Cancer Center
8 // Cedars Sinai Medical Center, Los Angeles, California 90048
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 
25 #include "Chemistry.hpp"
26 #include "Ion.hpp"
28 #include <cstring>
30 #include "boost/thread/thread.hpp"
31 #include "boost/thread/barrier.hpp"
32 
33 
34 using namespace pwiz::util;
35 using namespace pwiz::chemistry;
36 
37 
38 ostream* os_ = 0;
39 
40 
42 {
43  MassAbundance ma(1, 2);
44  MassAbundance ma2(1, 4);
45  unit_assert(ma != ma2);
46  ma2.abundance = 2;
47  unit_assert(ma == ma2);
48 }
49 
50 
52 {
53  const char* formula;
54  int numC, numH, numN, numO, numS, num13C, num2H, num15N, num18O;
55  double monoMass;
56  double avgMass;
57 };
58 
60 {
61  { "C1H2N3O4S5", 1, 2, 3, 4, 5, 0, 0, 0, 0, 279.864884, 280.36928 },
62  { "C1H 2N3O4 S5", 1, 2, 3, 4, 5, 0, 0, 0, 0, 279.864884, 280.36928 },
63  { "H-42", 0, -42, 0, 0, 0, 0, 0, 0, 0, -42.328651, -42.333512 },
64  { "N2C-1", -1, 0, 2, 0, 0, 0, 0, 0, 0, 28.006148 - 12, 28.013486 - 12.0107 },
65  { "C39H67N11O10", 39, 67, 11, 10, 0, 0, 0, 0, 0, 849.507238, 850.01698 },
66  { "C3H7N1O2Se1", 3, 7, 1, 2, 0, 0, 0, 0, 0, 168.9642, 168.0532 },
67  { "_13C1_2H2_15N3_18O4", 0, 0, 0, 0, 0, 1, 2, 3, 4, 134.0285, 134.0285 },
68  { "C-1 _13C1 _2H2 H-2 N-3 _15N3 _18O4 O-4", -1, -2, -3, -4, 0, 1, 2, 3, 4, 14.024, 13.984 },
69  { "C-1_13C1_2H2H-2N-3_15N3_18O4O-4", -1, -2, -3, -4, 0, 1, 2, 3, 4, 14.024, 13.984 },
70 };
71 
72 const int testFormulaDataSize = sizeof(testFormulaData)/sizeof(TestFormula);
73 
75 {
76  for (int i=0; i < testFormulaDataSize; ++i)
77  {
78  const TestFormula& testFormula = testFormulaData[i];
79  Formula formula(testFormula.formula);
80 
81  const double EPSILON = 0.001;
82 
83  if (os_) *os_ << formula << " " << formula.monoisotopicMass() << " " << formula.molecularWeight() << endl;
84  unit_assert_equal(formula.monoisotopicMass(), testFormula.monoMass, EPSILON);
85  unit_assert_equal(formula.molecularWeight(), testFormula.avgMass, EPSILON);
86 
87  formula[Element::C] += 2;
88  if (os_) *os_ << formula << " " << formula.monoisotopicMass() << " " << formula.molecularWeight() << endl;
89  unit_assert_equal(formula.monoisotopicMass(), testFormula.monoMass+24, EPSILON);
90  unit_assert_equal(formula.molecularWeight(), testFormula.avgMass+12.0107*2, EPSILON);
91 
92  //const Formula& constFormula = formula;
93  //constFormula[Element::C] = 1; // this won't compile, by design
94 
95  // test copy constructor
96  Formula formula2 = formula;
97  formula2[Element::C] -= 2;
98  if (os_) *os_ << "formula: " << formula << endl;
99  if (os_) *os_ << "formula2: " << formula2 << endl;
100  unit_assert_equal(formula.monoisotopicMass(), testFormula.monoMass+24, EPSILON);
101  unit_assert_equal(formula2.monoisotopicMass(), testFormula.monoMass, EPSILON);
102 
103  // test operator=
104  formula = formula2;
105  if (os_) *os_ << "formula: " << formula << endl;
106  if (os_) *os_ << "formula2: " << formula2 << endl;
107  unit_assert_equal(formula.monoisotopicMass(), formula2.monoisotopicMass(), EPSILON);
108 
109  // test operator==
110  unit_assert(formula == testFormula.formula); // implicit construction from string
111  unit_assert(formula == formula2);
112  formula2[Element::C] += 4; // test difference in CHONSP
113  unit_assert(formula != formula2);
114  formula2[Element::C] -= 4;
115  unit_assert(formula == formula2);
116  formula2[Element::U] += 2; // test difference outside CHONSP
117  unit_assert(formula != formula2);
118  formula2[Element::U] -= 2;
119  unit_assert(formula == formula2);
120 
121  // test data()
122  Formula::Map data = formula.data();
123  if (os_)
124  {
125  *os_ << "map: ";
126  for (Formula::Map::iterator it=data.begin(), end=data.end(); it!=end; ++it)
127  *os_ << it->first << it->second << " ";
128  *os_ << "\n";
129  }
130  }
131 }
132 
133 
135 {
136  using namespace Element;
137 
138  Formula water("H2O1");
139  Formula a("C1 H2 N3 O4 S5");
140 
141  a *= 2;
142  unit_assert(a[C]==2 && a[H]==4 && a[N]==6 && a[O]==8 && a[S]==10);
143  a += water;
144  unit_assert(a[H]==6 && a[O]==9);
145  a -= water;
146  unit_assert(a[C]==2 && a[H]==4 && a[N]==6 && a[O]==8 && a[S]==10);
147  a += 2*water;
148  unit_assert(a[H]==8 && a[O]==10);
149  a = (a - water*2);
150  unit_assert(a[C]==2 && a[H]==4 && a[N]==6 && a[O]==8 && a[S]==10);
151  a = water + water;
152  unit_assert(a[H]==4 && a[O]==2);
153  if (os_) *os_ << "water: " << a-water << endl;
154 }
155 
156 
157 void testInfo()
158 {
159  if (os_)
160  {
161  for (Element::Type e=Element::H; e<=Element::Ca; e=(Element::Type)(e+1))
162  *os_ << Element::Info::record(e).symbol << " " << Element::Info::record(e).atomicNumber << endl;
163  *os_ << endl;
164  }
165 
166  unit_assert(Element::Info::record(Element::C).atomicNumber == 6);
167  unit_assert(Element::Info::record("C").atomicNumber == 6);
168  unit_assert(Element::Info::record(Element::U).atomicNumber == 92);
169  unit_assert(Element::Info::record("U").atomicNumber == 92);
170  unit_assert(Element::Info::record(Element::Uuh).atomicNumber == 116);
171  unit_assert(Element::Info::record("Uuh").atomicNumber == 116);
172 
174  runtime_error,
175  "[chemistry::text2enum()] Error translating symbol foo");
176 }
177 
178 
180 {
181  if (os_)
182  {
183  *os_ << "Sulfur isotopes: " << Element::Info::record(Element::S).isotopes.size() << endl
184  << Element::Info::record(Element::S).isotopes;
185  }
186 }
187 
188 
190 {
191  Formula formula("Si6C12H36O6");
192 
193  if (os_)
194  {
195  *os_ << "polysiloxane:\n"
196  << formula << " "
197  << formula.monoisotopicMass() << " "
198  << formula.molecularWeight() << endl
199  << "ion: " << Ion::mz(formula.monoisotopicMass(), 1) << endl;
200  }
201 }
202 
203 
204 void testThreadSafetyWorker(boost::barrier* testBarrier, int& result)
205 {
206  testBarrier->wait(); // wait until all threads have started
207 
208  try
209  {
211  testFormula();
213  testInfo();
214  infoExample();
216  result = 0;
217  return;
218  }
219  catch (exception& e)
220  {
221  cerr << "Exception in worker thread: " << e.what() << endl;
222  }
223  catch (...)
224  {
225  cerr << "Unhandled exception in worker thread." << endl;
226  }
227  result = 1;
228 }
229 
230 void testThreadSafety(const int& testThreadCount)
231 {
232  boost::barrier testBarrier(testThreadCount);
233  boost::thread_group testThreadGroup;
234  vector<int> results(testThreadCount, 0);
235  for (int i=0; i < testThreadCount; ++i)
236  testThreadGroup.add_thread(new boost::thread(&testThreadSafetyWorker, &testBarrier, boost::ref(results[i])));
237  testThreadGroup.join_all();
238 
239  int failedThreads = std::accumulate(results.begin(), results.end(), 0);
240  if (failedThreads > 0)
241  throw runtime_error(lexical_cast<string>(failedThreads) + " thread failed");
242 }
243 
244 
245 int main(int argc, char* argv[])
246 {
247  TEST_PROLOG(argc, argv)
248 
249  if (argc>1 && !strcmp(argv[1],"-v")) os_ = &cout;
250  if (os_) *os_ << "ChemistryTest\n" << setprecision(12);
251 
252  try
253  {
254  //testThreadSafety(1); // does not test thread-safety of singleton initialization
255  testThreadSafety(2);
256  testThreadSafety(4);
257  testThreadSafety(8);
258  testThreadSafety(16);
259  }
260  catch (exception& e)
261  {
262  TEST_FAILED(e.what())
263  }
264  catch (...)
265  {
266  TEST_FAILED("Caught unknown exception.")
267  }
268 
270 }
271 
272 
void testMassAbundance()
#define unit_assert_throws_what(x, exception, whatStr)
Definition: unit.hpp:119
const TestFormula testFormulaData[]
class to represent a chemical formula
Definition: Chemistry.hpp:133
PWIZ_API_DECL const Record & record(Type type)
returns the amino acid&#39;s Record by type
void testFormula()
void testFormulaOperations()
#define TEST_EPILOG
Definition: unit.hpp:182
#define unit_assert_equal(x, y, epsilon)
Definition: unit.hpp:99
double molecularWeight() const
const int testFormulaDataSize
double monoisotopicMass() const
U
Definition: Chemistry.hpp:80
void testThreadSafetyWorker(boost::barrier *testBarrier, int &result)
int main(int argc, char *argv[])
ostream * os_
void infoExample()
N
Definition: Chemistry.hpp:80
Ca
Definition: Chemistry.hpp:80
double mz(double neutralMass, int protonDelta, int electronDelta=0, int neutronDelta=0)
Definition: Ion.hpp:78
void testPolysiloxane()
void testThreadSafety(const int &testThreadCount)
const char * formula
H
Definition: Chemistry.hpp:80
#define TEST_FAILED(x)
Definition: unit.hpp:176
S
Definition: Chemistry.hpp:80
#define TEST_PROLOG(argc, argv)
Definition: unit.hpp:174
pwiz::chemistry::Element::Type Type
O
Definition: Chemistry.hpp:80
C
Definition: Chemistry.hpp:80
struct for holding isotope information
Definition: Chemistry.hpp:51
#define unit_assert(x)
Definition: unit.hpp:85
std::map< Element::Type, int > Map
Definition: Chemistry.hpp:153
void testInfo()