Skills › Scientific Agent Skills › pyopenms › mass_calculator.py
mass_calculator.py
94 lines · skills/pyopenms/scripts/mass_calculator.py at 5f2bc03
1#!/usr/bin/env python3
2"""
3Mass & Chemistry Calculator
5Compute masses and isotope distributions for peptides (amino-acid sequences),
6empirical formulas, or both. Reports monoisotopic and average mass, m/z for a
7range of charge states, the molecular formula, and the theoretical isotope
8pattern.
10Usage:
11 python mass_calculator.py --peptide DFPIANGER
12 python mass_calculator.py --peptide "PEPTIDEM(Oxidation)K" --charges 1 2 3
13 python mass_calculator.py --formula C6H12O6 --isotopes 5
14 python mass_calculator.py --peptide DFPIANGER --isotopes 6 --csv iso.csv
15"""
17import argparse
18import sys
20try:
21 import pyopenms as ms
22except ImportError:
23 print("Error: pyopenms not installed. Install with: uv pip install pyopenms")
24 sys.exit(1)
26PROTON = 1.0072764665789
29def report_isotopes(formula, n, csv=None):
30 gen = ms.CoarseIsotopePatternGenerator(n)
31 dist = formula.getIsotopeDistribution(gen)
32 print(f"\nIsotope pattern (top {n}):")
33 rows = []
34 for iso in dist.getContainer():
35 mz = iso.getMZ()
36 prob = iso.getIntensity()
37 print(f" neutral mass {mz:.4f} Da rel.abundance {prob*100:6.2f}%")
38 rows.append((mz, prob))
39 if csv:
40 import csv as csvmod
41 with open(csv, "w", newline="") as fh:
42 w = csvmod.writer(fh)
43 w.writerow(["mass", "rel_abundance"])
44 w.writerows(rows)
45 print(f"Wrote {csv}")
48def main():
49 parser = argparse.ArgumentParser(description="Compute masses and isotope distributions.")
50 parser.add_argument("--peptide", help="Amino-acid sequence (TPP/OpenMS mod syntax allowed)")
51 parser.add_argument("--formula", help="Empirical formula, e.g. C6H12O6")
52 parser.add_argument("--charges", type=int, nargs="+", default=[1, 2, 3],
53 help="Charge states for m/z (default 1 2 3)")
54 parser.add_argument("--isotopes", type=int, default=0, help="Number of isotope peaks to report")
55 parser.add_argument("--negative", action="store_true", help="Report negative-mode m/z")
56 parser.add_argument("--csv", help="Write isotope pattern to CSV")
57 args = parser.parse_args()
58 if any(z <= 0 for z in args.charges):
59 parser.error("--charges must contain positive charge magnitudes")
61 if not args.peptide and not args.formula:
62 parser.error("provide --peptide and/or --formula")
64 sign = -1 if args.negative else 1
66 if args.peptide:
67 seq = ms.AASequence.fromString(args.peptide)
68 formula = seq.getFormula()
69 mono = seq.getMonoWeight()
70 avg = seq.getAverageWeight()
71 print(f"Peptide: {seq.toString()}")
72 print(f"Formula: {formula.toString()}")
73 print(f"Monoisotopic mass: {mono:.5f}")
74 print(f"Average mass: {avg:.5f}")
75 for z in args.charges:
76 mz = (mono + sign * z * PROTON) / z
77 print(f" [M{'+' if not args.negative else '-'}{z}H]{'+' if not args.negative else '-'} m/z = {mz:.5f} (z={z})")
78 if args.isotopes:
79 report_isotopes(formula, args.isotopes, args.csv if not args.formula else None)
81 if args.formula:
82 formula = ms.EmpiricalFormula(args.formula)
83 print(f"\nFormula: {formula.toString()}")
84 print(f"Monoisotopic mass: {formula.getMonoWeight():.5f}")
85 print(f"Average mass: {formula.getAverageWeight():.5f}")
86 for z in args.charges:
87 mz = (formula.getMonoWeight() + sign * z * PROTON) / z
88 print(f" m/z = {mz:.5f} (z={z})")
89 if args.isotopes:
90 report_isotopes(formula, args.isotopes, args.csv)
93if __name__ == "__main__":
94 main()