Skip to content

Commit 783ce48

Browse files
committed
getMass function for molecules, values mirror the function for formulas.
1 parent 1fbb515 commit 783ce48

2 files changed

Lines changed: 253 additions & 94 deletions

File tree

base/standard/src/main/java/org/openscience/cdk/tools/manipulator/AtomContainerManipulator.java

Lines changed: 182 additions & 77 deletions
Original file line numberDiff line numberDiff line change
@@ -206,102 +206,207 @@ public static double getTotalCharge(IAtomContainer atomContainer) {
206206
return charge;
207207
}
208208

209-
/**
210-
* Get the summed exact mass of all atoms in an AtomContainer. It
211-
* requires isotope information for all atoms to be set. Either set
212-
* this information using the {@link IsotopeFactory}, or use the
213-
* {@link MolecularFormulaManipulator#getMajorIsotopeMass(org.openscience.cdk.interfaces.IMolecularFormula)}
214-
* method, after converting the {@link IAtomContainer} to a
215-
* {@link IMolecularFormula} with the {@link MolecularFormulaManipulator}.
216-
*
217-
* @param atomContainer The IAtomContainer to manipulate
218-
* @return The summed exact mass of all atoms in this AtomContainer.
219-
* @see #getMolecularWeight(IAtomContainer)
220-
*/
221-
public static double getTotalExactMass(IAtomContainer atomContainer) {
222-
try {
209+
private static boolean hasIsotopeSpecified(IIsotope atom) {
210+
return atom.getMassNumber() != null && atom.getMassNumber() != 0;
211+
}
223212

224-
Isotopes isotopes = Isotopes.getInstance();
225-
double mass = 0.0;
226-
double hExactMass = isotopes.getMajorIsotope(1).getExactMass();
227-
for (IAtom atom : atomContainer.atoms()) {
228-
if (atom.getImplicitHydrogenCount() == null)
229-
throw new IllegalArgumentException("an atom had with unknown (null) implicit hydrogens");
230-
mass += atom.getExactMass();
231-
mass += atom.getImplicitHydrogenCount() * hExactMass;
232-
}
233-
return mass;
234-
} catch (IOException e) {
235-
throw new RuntimeException("Isotopes definitions could not be loaded", e);
213+
private static double getExactMass(IsotopeFactory isofact, IIsotope atom) {
214+
if (atom.getExactMass() != null)
215+
return atom.getExactMass();
216+
else if (atom.getMassNumber() != null)
217+
return isofact.getExactMass(getAtomicNum(atom),
218+
atom.getMassNumber());
219+
else
220+
return isofact.getMajorIsotopeMass(getAtomicNum(atom));
221+
}
222+
223+
private static double getMassOrAvg(IsotopeFactory isofact, IIsotope atom) {
224+
if (!hasIsotopeSpecified(atom))
225+
return isofact.getNaturalMass(atom);
226+
return getExactMass(isofact, atom);
227+
}
228+
229+
public static final int MolWeight = 0x1;
230+
public static final int AverageWeight = 0x2;
231+
public static final int MonoIsotopic = 0x3;
232+
public static final int MostAbundant = 0x4;
233+
234+
public static final Comparator<IIsotope> NAT_ABUN_COMP = new Comparator<IIsotope>() {
235+
@Override
236+
public int compare(IIsotope o1, IIsotope o2) {
237+
return -Double.compare(o1.getNaturalAbundance(),
238+
o2.getNaturalAbundance());
236239
}
240+
};
241+
242+
private static double getDistMass(IsotopeFactory isofact,
243+
IIsotope[] isos, int idx, int count) {
244+
if (count == 0)
245+
return 0;
246+
double frac = 100d;
247+
double res = 0;
248+
for (int i = 0; i < idx; i++)
249+
frac -= isos[i].getNaturalAbundance();
250+
double p = isos[idx].getNaturalAbundance() / frac;
251+
if (p >= 1.0)
252+
return count * isos[idx].getExactMass();
253+
double kMin = (count + 1) * (1 - p) - 1;
254+
double kMax = (count + 1) * (1 - p);
255+
if ((int) Math.ceil(kMin) == (int) Math.floor(kMax)) {
256+
int k = (int) kMax;
257+
res = (count - k) * getExactMass(isofact,
258+
isos[idx]);
259+
res += getDistMass(isofact, isos, idx + 1, k);
260+
}
261+
return res;
262+
}
263+
264+
private static int getImplHCount(IAtom atom) {
265+
Integer implh = atom.getImplicitHydrogenCount();
266+
if (implh == null)
267+
throw new IllegalArgumentException("An atom had 'null' implicit hydrogens!");
268+
return implh;
269+
}
270+
271+
private static int getAtomicNum(IElement atom) {
272+
Integer atno = atom.getAtomicNumber();
273+
if (atno == null)
274+
throw new IllegalArgumentException("An atom had 'null' atomic number!");
275+
return atno;
237276
}
238277

239278
/**
240-
* Returns the molecular mass of the IAtomContainer. For the calculation it
241-
* uses the masses of the isotope mixture using natural abundances.
279+
* Calculate the mass of a molecule, this function takes an optional
280+
* 'mass flavour' that switches the computation type, the flavours are:
281+
* <br>
282+
* <ul>
283+
* <li>{@link #MolWeight} (default) - use and isotopes the natural mass
284+
* unless a specific isotope is specified</li>
285+
* <li>{@link #AverageWeight} - use and isotopes the natural mass even
286+
* if a specific isotope is specified</li>
287+
* <li>{@link #MonoIsotopic} - use and isotopes the major isotope mass
288+
* even if a specific isotope is specified</li>
289+
* <li>{@link #MostAbundant} - use the distribution of isotopes
290+
* based on their abundance and select the most abundant. For example
291+
* C<sub>6</sub>Br<sub>6</sub> would have three <sup>79</sup>Br and
292+
* <sup>81</sup>Br because their abundance is 51 and 49%.
293+
* </ul>
242294
*
243-
* @param atomContainer
244-
* @cdk.keyword mass, molecular
245-
* @see #getMolecularWeight(IAtomContainer)
295+
* @param mol molecule to compute mass for
296+
* @param flav the mass flavour
297+
* @return the mass of the molecule
298+
* @see #MolWeight
299+
* @see #AverageWeight
300+
* @see #MonoIsotopic
301+
* @see #MostAbundant
246302
*/
247-
public static double getNaturalExactMass(IAtomContainer atomContainer) {
248-
try {
249-
Isotopes isotopes = Isotopes.getInstance();
250-
double hydgrogenMass = isotopes.getNaturalMass(Elements.HYDROGEN);
303+
public static double getMass(IAtomContainer mol, int flav) {
251304

252-
double mass = 0.0;
253-
for (final IAtom atom : atomContainer.atoms()) {
305+
final Isotopes isofact;
306+
try {
307+
isofact = Isotopes.getInstance();
308+
} catch (IOException e) {
309+
throw new IllegalStateException("Could not load Isotopes!");
310+
}
254311

255-
if (atom.getAtomicNumber() == null)
256-
throw new IllegalArgumentException("an atom had with unknown (null) atomic number");
257-
if (atom.getImplicitHydrogenCount() == null)
258-
throw new IllegalArgumentException("an atom had with unknown (null) implicit hydrogens");
312+
double mass = 0;
313+
int hcnt = 0;
259314

260-
mass += isotopes.getNaturalMass(Elements.ofNumber(atom.getAtomicNumber()).toIElement());
261-
mass += hydgrogenMass * atom.getImplicitHydrogenCount();
262-
}
263-
return mass;
315+
switch (flav & 0xf) {
316+
case MolWeight:
317+
for (IAtom atom : mol.atoms()) {
318+
mass += getMassOrAvg(isofact, atom);
319+
hcnt += getImplHCount(atom);
320+
}
321+
mass += hcnt * isofact.getNaturalMass(1);
322+
break;
323+
case AverageWeight:
324+
for (IAtom atom : mol.atoms()) {
325+
mass += isofact.getNaturalMass(getAtomicNum(atom));
326+
hcnt += getImplHCount(atom);
327+
}
328+
mass += hcnt * isofact.getNaturalMass(1);
329+
break;
330+
case MonoIsotopic:
331+
for (IAtom atom : mol.atoms()) {
332+
mass += getExactMass(isofact, atom);
333+
hcnt += getImplHCount(atom);
334+
}
335+
mass += hcnt * isofact.getMajorIsotopeMass(1);
336+
break;
337+
case MostAbundant:
338+
int[] mf = new int[128];
339+
for (IAtom atom : mol.atoms()) {
340+
if (hasIsotopeSpecified(atom))
341+
mass += getExactMass(isofact, atom);
342+
else
343+
mf[getAtomicNum(atom)]++;
344+
mf[1] += atom.getImplicitHydrogenCount();
345+
}
264346

265-
} catch (IOException e) {
266-
throw new RuntimeException("Isotopes definitions could not be loaded", e);
347+
for (int atno = 0; atno < mf.length; atno++) {
348+
if (mf[atno] == 0)
349+
continue;
350+
IIsotope[] isotopes = isofact.getIsotopes(atno);
351+
Arrays.sort(isotopes, NAT_ABUN_COMP);
352+
mass += getDistMass(isofact, isotopes, 0, mf[atno]);
353+
}
354+
break;
267355
}
356+
return mass;
268357
}
269358

270359
/**
271-
* Calculate the molecular weight of a molecule.
360+
* Calculate the mass of a molecule, this function takes an optional
361+
* 'mass flavour' that switches the computation type, the flavours are:
362+
* <br>
363+
* <ul>
364+
* <li>{@link #MolWeight} (default) - use and isotopes the natural mass
365+
* unless a specific isotope is specified</li>
366+
* <li>{@link #AverageWeight} - use and isotopes the natural mass even
367+
* if a specific isotope is specified</li>
368+
* <li>{@link #MonoIsotopic} - use and isotopes the major isotope mass
369+
* even if a specific isotope is specified</li>
370+
* <li>{@link #MostAbundant} - use the distribution of isotopes
371+
* based on their abundance and select the most abundant. For example
372+
* C<sub>6</sub>Br<sub>6</sub> would have three <sup>79</sup>Br and
373+
* <sup>81</sup>Br because their abundance is 51 and 49%.
374+
* </ul>
272375
*
273-
* @param mol the molecule
274-
* @return the molecular weight
376+
* @param mol molecule to compute mass for
377+
* @return the mass of the molecule
378+
* @see #getMass(IAtomContainer, int)
379+
* @see #MolWeight
380+
* @see #AverageWeight
381+
* @see #MonoIsotopic
382+
* @see #MostAbundant
275383
*/
276-
public static double getMolecularWeight(IAtomContainer mol) {
277-
try {
278-
Isotopes isotopes = Isotopes.getInstance();
279-
double hmass = isotopes.getNaturalMass(Elements.HYDROGEN);
280-
double mw = 0.0;
281-
for (final IAtom atom : mol.atoms()) {
282-
if (atom.getAtomicNumber() == null || atom.getAtomicNumber() == 0)
283-
throw new IllegalArgumentException("An atom had with unknown (null) atomic number");
284-
if (atom.getImplicitHydrogenCount() == null)
285-
throw new IllegalArgumentException("An atom had with unknown (null) implicit hydrogens");
286-
mw += hmass * atom.getImplicitHydrogenCount();
287-
if (atom.getMassNumber() == null)
288-
mw += isotopes.getNaturalMass(atom);
289-
else if (atom.getExactMass() != null)
290-
mw += atom.getExactMass();
291-
else {
292-
IIsotope isotope = isotopes.getIsotope(atom.getSymbol(), atom.getMassNumber());
293-
if (isotope == null)
294-
mw += isotopes.getNaturalMass(atom);
295-
else
296-
mw += isotope.getExactMass();
297-
}
384+
public static double getMass(IAtomContainer mol) {
385+
return getMass(mol, MolWeight);
386+
}
298387

299-
}
300-
return mw;
388+
/**
389+
* @deprecated use {@link #getMass(IAtomContainer, int)} and
390+
* {@link #MonoIsotopic}
391+
*/
392+
public static double getTotalExactMass(IAtomContainer mol) {
393+
return getMass(mol, MonoIsotopic);
394+
}
301395

302-
} catch (IOException e) {
303-
throw new RuntimeException("Isotopes definitions could not be loaded", e);
304-
}
396+
/**
397+
* @deprecated use {@link #getMass(IAtomContainer, int)} and
398+
* {@link #AverageWeight}
399+
*/
400+
public static double getNaturalExactMass(IAtomContainer mol) {
401+
return getMass(mol, AverageWeight);
402+
}
403+
404+
/**
405+
* @deprecated use {@link #getMass(IAtomContainer, int)} and
406+
* {@link #MolWeight}
407+
*/
408+
public static double getMolecularWeight(IAtomContainer mol) {
409+
return getMass(mol, MolWeight);
305410
}
306411

307412
/**

base/test-standard/src/test/java/org/openscience/cdk/tools/manipulator/AtomContainerManipulatorTest.java

Lines changed: 71 additions & 17 deletions
Original file line numberDiff line numberDiff line change
@@ -18,23 +18,6 @@
1818
*/
1919
package org.openscience.cdk.tools.manipulator;
2020

21-
import static org.hamcrest.CoreMatchers.is;
22-
import static org.hamcrest.CoreMatchers.not;
23-
import static org.hamcrest.CoreMatchers.sameInstance;
24-
import static org.hamcrest.number.IsCloseTo.closeTo;
25-
import static org.junit.Assert.assertFalse;
26-
import static org.junit.Assert.assertNotNull;
27-
import static org.junit.Assert.assertNull;
28-
import static org.junit.Assert.assertThat;
29-
import static org.junit.Assert.assertTrue;
30-
31-
import java.io.IOException;
32-
import java.io.InputStream;
33-
import java.util.Collection;
34-
import java.util.HashMap;
35-
import java.util.List;
36-
import java.util.Map;
37-
3821
import org.hamcrest.CoreMatchers;
3922
import org.junit.Assert;
4023
import org.junit.Before;
@@ -70,6 +53,18 @@
7053
import org.openscience.cdk.templates.TestMoleculeFactory;
7154
import org.openscience.cdk.tools.CDKHydrogenAdder;
7255

56+
import java.io.IOException;
57+
import java.io.InputStream;
58+
import java.util.Collection;
59+
import java.util.HashMap;
60+
import java.util.List;
61+
import java.util.Map;
62+
63+
import static org.hamcrest.CoreMatchers.*;
64+
import static org.hamcrest.number.IsCloseTo.closeTo;
65+
import static org.junit.Assert.*;
66+
import static org.openscience.cdk.tools.manipulator.AtomContainerManipulator.*;
67+
7368
/**
7469
* @cdk.module test-standard
7570
*/
@@ -1339,4 +1334,63 @@ static void assertRemoveH(String smiIn, String smiExp) throws Exception {
13391334

13401335
assertThat(smiAct, is(smiExp));
13411336
}
1337+
1338+
@Test public void getMassC6Br6() throws InvalidSmilesException {
1339+
IChemObjectBuilder bldr = SilentChemObjectBuilder.getInstance();
1340+
SmilesParser smipar = new SmilesParser(bldr);
1341+
IAtomContainer mol = smipar.parseSmiles("Brc1c(Br)c(Br)c(Br)c(Br)c1Br");
1342+
assertThat(AtomContainerManipulator.getMass(mol, MolWeight),
1343+
closeTo(551.485, 0.001));
1344+
assertThat(AtomContainerManipulator.getMass(mol, AverageWeight),
1345+
closeTo(551.485, 0.001));
1346+
assertThat(AtomContainerManipulator.getMass(mol, MonoIsotopic),
1347+
closeTo(545.510, 0.001));
1348+
assertThat(AtomContainerManipulator.getMass(mol, MostAbundant),
1349+
closeTo(551.503, 0.001));
1350+
}
1351+
1352+
@Test public void getMassCranbin() {
1353+
IChemObjectBuilder bldr = SilentChemObjectBuilder.getInstance();
1354+
IAtomContainer mol =
1355+
MolecularFormulaManipulator.getAtomContainer("C202H315N55O64S6",
1356+
bldr);
1357+
assertThat(AtomContainerManipulator.getMass(mol, MolWeight),
1358+
closeTo(4730.397, 0.001));
1359+
assertThat(AtomContainerManipulator.getMass(mol, AverageWeight),
1360+
closeTo(4730.397, 0.001));
1361+
assertThat(AtomContainerManipulator.getMass(mol, MonoIsotopic),
1362+
closeTo(4727.140, 0.001));
1363+
assertThat(AtomContainerManipulator.getMass(mol, MostAbundant),
1364+
closeTo(4729.147, 0.001));
1365+
}
1366+
1367+
@Test public void getMassCranbinSpecIsotopes() {
1368+
IChemObjectBuilder bldr = SilentChemObjectBuilder.getInstance();
1369+
IAtomContainer mol =
1370+
MolecularFormulaManipulator.getAtomContainer("[12]C200[13]C2[1]H315[14]N55[16]O64[32]S6",
1371+
bldr);
1372+
assertThat(AtomContainerManipulator.getMass(mol, MolWeight),
1373+
closeTo(4729.147, 0.001));
1374+
assertThat(AtomContainerManipulator.getMass(mol, AverageWeight),
1375+
closeTo(4730.397, 0.001));
1376+
assertThat(AtomContainerManipulator.getMass(mol, MonoIsotopic),
1377+
closeTo(4729.147, 0.001));
1378+
assertThat(AtomContainerManipulator.getMass(mol, MostAbundant),
1379+
closeTo(4729.147, 0.001));
1380+
}
1381+
1382+
@Test public void getMassCranbinMixedSpecIsotopes() {
1383+
IChemObjectBuilder bldr = SilentChemObjectBuilder.getInstance();
1384+
IAtomContainer mol =
1385+
MolecularFormulaManipulator.getAtomContainer("C200[13]C2H315N55O64S6",
1386+
bldr);
1387+
assertThat(AtomContainerManipulator.getMass(mol, MolWeight),
1388+
closeTo(4732.382, 0.001));
1389+
assertThat(AtomContainerManipulator.getMass(mol, AverageWeight),
1390+
closeTo(4730.397, 0.001));
1391+
assertThat(AtomContainerManipulator.getMass(mol, MonoIsotopic),
1392+
closeTo(4729.147, 0.001));
1393+
assertThat(AtomContainerManipulator.getMass(mol, MostAbundant),
1394+
closeTo(4731.154, 0.001));
1395+
}
13421396
}

0 commit comments

Comments
 (0)