ProteoWizard
IsotopeTableTest.cpp
Go to the documentation of this file.
1//
2// $Id$
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
24#include "IsotopeTable.hpp"
27#include <cstring>
28
29
30using namespace pwiz::util;
31using namespace pwiz::chemistry;
32
33
34ostream* os_ = 0;
35
36
38{
39 return a.abundance > b.abundance;
40}
41
42
43bool hasLessMass(const MassAbundance& a, const MassAbundance& b)
44{
45 return a.mass < b.mass;
46}
47
48
49void test1()
50{
52 md.push_back(MassAbundance(10, 1));
53
54 IsotopeTable table(md, 10, 0);
55
56 for (int i=1; i<=10; i++)
57 {
58 MassDistribution temp = table.distribution(i);
59 unit_assert(temp.size() == 1);
60 unit_assert(temp[0] == MassAbundance(i*10, 1));
61 //cout << i << " atoms:\n" << temp << endl;
62 }
63}
64
65
66void test2()
67{
68 const double p0 = .9;
69 const double p1 = 1 - p0;
70
72 md.push_back(MassAbundance(10, p0));
73 md.push_back(MassAbundance(11, p1));
74
75 IsotopeTable table(md, 10, 0);
76
77/*
78 for (int i=0; i<=10; i++)
79 cout << i << " atoms:\n" << table.distribution(i) << endl;
80*/
81
82 // test manually for 1 atom
83
85 unit_assert(test1.size() == 2);
86 unit_assert(test1[0] == md[0]);
87 unit_assert(test1[1] == md[1]);
88
89 // test manually for 10 atoms
90
91 const int n = 10;
92 MassDistribution good10;
93 double abundance = pow(p0, n);
94 double mass = 100;
95
96 for (int k=0; k<=n; k++)
97 {
98 good10.push_back(MassAbundance(mass, abundance));
99 abundance *= p1/p0*(n-k)/(k+1);
100 mass += 1;
101 }
102
103 sort(good10.begin(), good10.end(), hasGreaterAbundance);
104
105 MassDistribution test10 = table.distribution(10);
106 sort(test10.begin(), test10.end(), hasGreaterAbundance);
107
108 unit_assert((int)test10.size() == n+1);
109
110 for (int k=0; k<=n; k++)
111 unit_assert_equal(test10[k].abundance, good10[k].abundance, 1e-15);
112
113 // test cutoff
114
115 IsotopeTable table_cutoff(md, 10, 1e-8);
116 unit_assert(table_cutoff.distribution(10).size() == 9);
117}
118
119
121{
122 unit_assert(test.size() == good.size());
123 for (unsigned int i=0; i<test.size(); i++)
124 {
125 unit_assert_equal(test[i].mass, good[i].mass, 1e-12);
126 unit_assert_equal(test[i].abundance, good[i].abundance, 1e-12);
127 }
128}
129
130
131void test3()
132{
133 const double p0 = .9;
134 const double p1 = .09;
135 const double p2 = 1 - (p0 + p1);
136
137 const double m0 = 10;
138 const double m1 = 11;
139 const double m2 = 12.33;
140
142 md.push_back(MassAbundance(m0, p0));
143 md.push_back(MassAbundance(m1, p1));
144 md.push_back(MassAbundance(m2, p2));
145
146// cout << "test3 distribution:\n" << md << endl;
147
148 IsotopeTable table(md, 10, 1e-5);
149
150 // compare distribution for 1 atom
151 compare(table.distribution(1), md);
152
153 // compare distribution for 2 atoms
154 MassDistribution good3_2;
155 good3_2.push_back(MassAbundance(m0*2, p0*p0));
156 good3_2.push_back(MassAbundance(m0+m1, p0*p1*2));
157 good3_2.push_back(MassAbundance(m0+m2, p0*p2*2));
158 good3_2.push_back(MassAbundance(m1+m2, p1*p2*2));
159 good3_2.push_back(MassAbundance(m1+m1, p1*p1));
160 good3_2.push_back(MassAbundance(m2+m2, p2*p2));
161 sort(good3_2.begin(), good3_2.end(), hasGreaterAbundance);
162
163 MassDistribution test3_2 = table.distribution(2);
164 sort(test3_2.begin(), test3_2.end(), hasGreaterAbundance);
165
166// cout << "good:\n" << good3_2 << endl;
167// cout << "test:\n" << test3_2 << endl;
168
169 compare(test3_2, good3_2);
170}
171
172
173void test4()
174{
175 const double p0 = .9;
176 const double p1 = .09;
177 const double p2 = .009;
178 const double p3 = 1 - (p0 + p1 + p2);
179
180 const double m0 = 10;
181 const double m1 = 11;
182 const double m2 = 12.33;
183 const double m3 = 13.13;
184
186 md.push_back(MassAbundance(m0, p0));
187 md.push_back(MassAbundance(m1, p1));
188 md.push_back(MassAbundance(m2, p2));
189 md.push_back(MassAbundance(m3, p3));
190
191 cout << "test4 distribution:\n" << md << endl;
192
193 IsotopeTable table(md, 10, 1e-5);
194
195 compare(md, table.distribution(1));
196
197 MassDistribution test4_2 = table.distribution(2);
198
199 cout << "2 atoms:\n" << test4_2 << endl;
200}
201
202
204{
205 IsotopeTable table(Element::Info::record(Element::Se).isotopes, 10, 1e-10);
206
207 cout << table << endl;
208
209 MassDistribution dist10 = table.distribution(10);
210 cout << "distribution: " << dist10 << endl;
211}
212
213
214int main(int argc, char* argv[])
215{
216 TEST_PROLOG(argc, argv)
217
218 try
219 {
220 if (argc>1 && !strcmp(argv[1],"-v")) os_ = &cout;
221 if (os_) *os_ << "IsotopeTableTest\n" << setprecision(12);
222 test1();
223 test2();
224 test3();
225 //test4();
226 //testSulfur();
227 }
228 catch (exception& e)
229 {
230 TEST_FAILED(e.what())
231 }
232 catch (...)
233 {
234 TEST_FAILED("Caught unknown exception.")
235 }
236
238}
239
void test2()
int main(int argc, char *argv[])
void compare(const MassDistribution &test, const MassDistribution &good)
void test1()
void test4()
void testSulfur()
void test3()
bool hasLessMass(const MassAbundance &a, const MassAbundance &b)
bool hasGreaterAbundance(const MassAbundance &a, const MassAbundance &b)
ostream * os_
Kernel k
Class representing a table of isotope distributions for collections of multiple atoms of a single ele...
MassDistribution distribution(int atomCount) const
PWIZ_API_DECL const Record & record(Type type)
retrieve the record for an element
std::vector< MassAbundance > MassDistribution
struct for holding isotope distribution
Definition Chemistry.hpp:66
struct for holding isotope information
Definition Chemistry.hpp:52
#define unit_assert(x)
Definition unit.hpp:85
#define TEST_EPILOG
Definition unit.hpp:183
#define TEST_FAILED(x)
Definition unit.hpp:177
#define unit_assert_equal(x, y, epsilon)
Definition unit.hpp:99
#define TEST_PROLOG(argc, argv)
Definition unit.hpp:175