1use chematic_core::{Atom, AtomIdx, BondOrder, Element, Molecule, MoleculeBuilder};
13
14use crate::error::MolParseError;
15
16const MAX_ATOMS: usize = 100_000;
18
19const MAX_BONDS: usize = 200_000;
21
22#[derive(Debug, Clone, PartialEq, Eq, Default)]
24pub struct MolMetadata {
25 pub name: String,
27 pub comment: String,
29}
30
31impl MolMetadata {
32 pub fn with_name(mut self, name: &str) -> Self {
34 self.name = name.to_owned();
35 self
36 }
37
38 pub fn with_comment(mut self, comment: &str) -> Self {
40 self.comment = comment.to_owned();
41 self
42 }
43}
44
45fn decode_charge(code: i8) -> i8 {
51 match code {
52 1 => 3,
53 2 => 2,
54 3 => 1,
55 4 => 0, 5 => -1,
57 6 => -2,
58 7 => -3,
59 _ => 0,
60 }
61}
62
63fn encode_charge(charge: i8) -> u8 {
65 match charge {
66 3 => 1,
67 2 => 2,
68 1 => 3,
69 -1 => 5,
70 -2 => 6,
71 -3 => 7,
72 _ => 0,
73 }
74}
75
76fn parse_field3(
85 line: &str,
86 start: usize,
87 line_num: usize,
88 make_err: impl Fn(usize, String) -> MolParseError,
89) -> Result<usize, MolParseError> {
90 let field = line
91 .get(start..start + 3)
92 .ok_or_else(|| make_err(line_num, format!("line too short at column {start}")))?;
93 field
94 .trim()
95 .parse::<usize>()
96 .map_err(|_| make_err(line_num, format!("cannot parse integer from '{field}'")))
97}
98
99#[allow(clippy::type_complexity)]
105pub fn parse_mol_with_coords(
106 input: &str,
107) -> Result<(Molecule, MolMetadata, Vec<(f64, f64)>), MolParseError> {
108 let mut lines = input.lines().enumerate().map(|(i, l)| (i + 1, l));
110 let mut next_line = || lines.next().ok_or(MolParseError::UnexpectedEnd);
111
112 let name = next_line()?.1.to_string();
115 next_line()?; let comment = next_line()?.1.to_string();
117
118 let metadata = MolMetadata { name, comment };
119
120 let (counts_lineno, counts_line) = next_line()?;
123
124 if !counts_line.contains("V2000") {
126 return Err(MolParseError::InvalidCountLine {
127 line: counts_lineno,
128 detail: "missing V2000 version tag".to_string(),
129 });
130 }
131
132 let make_count_err = |ln: usize, d: String| MolParseError::InvalidCountLine {
133 line: ln,
134 detail: d,
135 };
136
137 let natoms = parse_field3(counts_line, 0, counts_lineno, make_count_err)?;
138 let nbonds = parse_field3(counts_line, 3, counts_lineno, make_count_err)?;
139
140 if natoms > MAX_ATOMS {
141 return Err(MolParseError::InvalidCountLine {
142 line: counts_lineno,
143 detail: format!(
144 "atom count {} exceeds maximum allowed {}",
145 natoms, MAX_ATOMS
146 ),
147 });
148 }
149
150 if nbonds > MAX_BONDS {
151 return Err(MolParseError::InvalidCountLine {
152 line: counts_lineno,
153 detail: format!(
154 "bond count {} exceeds maximum allowed {}",
155 nbonds, MAX_BONDS
156 ),
157 });
158 }
159
160 let mut builder = MoleculeBuilder::new();
163 let mut coords: Vec<(f64, f64)> = Vec::with_capacity(natoms);
164 let make_atom_err = |ln: usize, d: String| MolParseError::InvalidAtomLine {
165 line: ln,
166 detail: d,
167 };
168
169 for atom_i in 0..natoms {
170 let (raw_lineno, atom_line) = next_line()?;
171
172 let x: f64 = atom_line
174 .get(0..10)
175 .and_then(|s| s.trim().parse().ok())
176 .unwrap_or(0.0);
177 let y: f64 = atom_line
178 .get(10..20)
179 .and_then(|s| s.trim().parse().ok())
180 .unwrap_or(0.0);
181 coords.push((x, y));
182
183 let sym = atom_line
186 .get(31..34)
187 .ok_or_else(|| {
188 make_atom_err(
189 raw_lineno,
190 format!("atom line {atom_i} too short for element field"),
191 )
192 })?
193 .trim();
194
195 let element = Element::from_symbol(sym).ok_or_else(|| MolParseError::UnknownElement {
196 symbol: sym.to_string(),
197 line: raw_lineno,
198 })?;
199
200 let charge = atom_line
202 .get(36..39)
203 .map(|ccc| decode_charge(ccc.trim().parse().unwrap_or(0)))
204 .unwrap_or(0);
205
206 let mut atom = Atom::new(element);
207 atom.charge = charge;
208 builder.add_atom(atom);
209 }
210
211 let make_bond_err = |ln: usize, d: String| MolParseError::InvalidBondLine {
214 line: ln,
215 detail: d,
216 };
217
218 for bond_i in 0..nbonds {
219 let (raw_lineno, bond_line) = next_line()?;
220
221 let a1_raw = parse_field3(bond_line, 0, raw_lineno, make_bond_err)?;
222 let a2_raw = parse_field3(bond_line, 3, raw_lineno, make_bond_err)?;
223 let btype_raw = parse_field3(bond_line, 6, raw_lineno, make_bond_err)?;
224
225 if a1_raw == 0 || a2_raw == 0 {
226 return Err(MolParseError::InvalidBondLine {
227 line: raw_lineno,
228 detail: format!("bond {bond_i}: atom indices are 1-based; got {a1_raw}/{a2_raw}"),
229 });
230 }
231
232 let a1 = AtomIdx((a1_raw - 1) as u32);
233 let a2 = AtomIdx((a2_raw - 1) as u32);
234
235 let stereo_raw: usize = if bond_line.len() >= 12 {
237 parse_field3(bond_line, 9, raw_lineno, make_bond_err).unwrap_or(0)
238 } else {
239 0
240 };
241
242 let order = match btype_raw {
243 1 => match stereo_raw {
244 1 | 4 => BondOrder::Up,
245 6 => BondOrder::Down,
246 _ => BondOrder::Single,
247 },
248 2 => BondOrder::Double,
249 3 => BondOrder::Triple,
250 4 => BondOrder::Aromatic,
251 5 => BondOrder::QuerySingleOrDouble,
252 6 => BondOrder::QuerySingleOrAromatic,
253 7 => BondOrder::QueryDoubleOrAromatic,
254 8 => BondOrder::QueryAny,
255 _ => BondOrder::Single,
256 };
257
258 builder
259 .add_bond(a1, a2, order)
260 .map_err(|e| MolParseError::InvalidBondLine {
261 line: raw_lineno,
262 detail: format!("bond {bond_i}: {e}"),
263 })?;
264 }
265
266 for (_, l) in lines.by_ref() {
268 if l.trim_start().starts_with("M END") {
269 break;
270 }
271 }
272
273 Ok((builder.build(), metadata, coords))
274}
275
276pub fn parse_mol(input: &str) -> Result<(Molecule, MolMetadata), MolParseError> {
281 parse_mol_with_coords(input).map(|(mol, meta, _coords)| (mol, meta))
282}
283
284#[allow(clippy::type_complexity)]
291pub fn parse_sdf_with_coords(
292 input: &str,
293) -> Result<Vec<(Molecule, MolMetadata, Vec<(f64, f64)>)>, MolParseError> {
294 let mut result = Vec::new();
297 let mut remaining = input;
298 loop {
299 while let Some(rest) = remaining
301 .strip_prefix("\r\n")
302 .or_else(|| remaining.strip_prefix('\n'))
303 {
304 remaining = rest;
305 }
306 if remaining.is_empty() {
307 break;
308 }
309
310 let mut byte_offset = 0usize;
312 let (end_byte, after_delim) = loop {
313 let rest = &remaining[byte_offset..];
314 match rest.find('\n') {
315 Some(nl) => {
316 let line = rest[..nl].trim_end_matches('\r');
317 if line == "$$$$" {
318 break (byte_offset, &remaining[byte_offset + nl + 1..]);
319 }
320 byte_offset += nl + 1;
321 }
322 None => {
323 if rest.trim_end_matches('\r') == "$$$$" {
324 break (byte_offset, "");
325 }
326 break (remaining.len(), "");
327 }
328 }
329 };
330
331 let block = &remaining[..end_byte];
332 remaining = after_delim;
333 if block.trim().is_empty() {
334 continue;
335 }
336
337 let (mol, meta, coords) = parse_mol_with_coords(block)?;
338 result.push((mol, meta, coords));
339 }
340 Ok(result)
341}
342
343pub fn write_mol(mol: &Molecule, metadata: &MolMetadata) -> String {
353 write_mol_with_coords(mol, metadata, &[])
354}
355
356pub fn write_mol_with_coords(
361 mol: &Molecule,
362 metadata: &MolMetadata,
363 coords: &[(f64, f64)],
364) -> String {
365 let mut out = String::new();
366
367 out.push_str(&metadata.name);
369 out.push('\n');
370 out.push_str(" chematic\n");
371 out.push_str(&metadata.comment);
372 out.push('\n');
373
374 let natoms = mol.atom_count();
376 let nbonds = mol.bond_count();
377 out.push_str(&format!(
378 "{:>3}{:>3} 0 0 0 0 0 0 0 0999 V2000\n",
379 natoms, nbonds
380 ));
381
382 for (idx, atom) in mol.atoms() {
384 let sym = atom.element.symbol();
385 let charge_code = encode_charge(atom.charge);
386 let (x, y) = coords.get(idx.0 as usize).copied().unwrap_or((0.0, 0.0));
387 out.push_str(&format!(
388 "{:>10.4}{:>10.4}{:>10.4} {:<3} 0{:>3} 0 0 0 0 0 0 0 0 0\n",
389 x, y, 0.0_f64, sym, charge_code,
390 ));
391 }
392
393 for (_idx, bond) in mol.bonds() {
395 let a1 = bond.atom1.0 + 1; let a2 = bond.atom2.0 + 1;
397 let btype = match bond.order {
398 BondOrder::Single | BondOrder::Up | BondOrder::Down | BondOrder::Dative => 1,
399 BondOrder::Double => 2,
400 BondOrder::Triple => 3,
401 BondOrder::Aromatic => 4,
402 BondOrder::QuerySingleOrDouble => 5,
403 BondOrder::QuerySingleOrAromatic => 6,
404 BondOrder::QueryDoubleOrAromatic => 7,
405 BondOrder::QueryAny | BondOrder::Zero => 8,
406 BondOrder::Quadruple => 4,
407 };
408 out.push_str(&format!("{:>3}{:>3}{:>3} 0\n", a1, a2, btype));
409 }
410
411 out.push_str("M END\n");
413
414 out
415}
416
417#[allow(clippy::type_complexity)]
427pub fn write_sdf(records: &[(&Molecule, &MolMetadata, &[(f64, f64)])]) -> String {
428 let mut out = String::new();
429 for (mol, meta, coords) in records {
430 out.push_str(&write_mol_with_coords(mol, meta, coords));
431 out.push_str("$$$$\n");
432 }
433 out
434}
435
436#[allow(clippy::type_complexity)]
451pub fn write_sdf_with_charges(
452 records: &[(&Molecule, &MolMetadata, &[(f64, f64)], &[f64])],
453) -> String {
454 let mut out = String::new();
455 for (mol, meta, coords, charges) in records {
456 out.push_str(&write_mol_with_coords(mol, meta, coords));
457 if !charges.is_empty() {
458 out.push_str("> <PARTIAL_CHARGES>\n");
459 let vals: Vec<String> = charges.iter().map(|q| format!("{q:.4}")).collect();
460 out.push_str(&vals.join(" "));
461 out.push_str("\n\n");
462 }
463 out.push_str("$$$$\n");
464 }
465 out
466}
467
468#[cfg(test)]
473mod tests {
474 use super::*;
475
476 const ETHANOL_MOL: &str = "\
478ethanol
479 chematic
480
481 3 2 0 0 0 0 0 0 0 0 0 V2000
482 0.0000 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0
483 1.5000 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0
484 3.0000 0.0000 0.0000 O 0 0 0 0 0 0 0 0 0 0 0 0
485 1 2 1 0
486 2 3 1 0
487M END
488";
489
490 #[test]
491 fn test_parse_ethanol_counts() {
492 let (mol, meta) = parse_mol(ETHANOL_MOL).expect("parse should succeed");
493 assert_eq!(mol.atom_count(), 3);
494 assert_eq!(mol.bond_count(), 2);
495 assert_eq!(meta.name, "ethanol");
496 }
497
498 #[test]
499 fn test_parse_elements() {
500 let (mol, _) = parse_mol(ETHANOL_MOL).expect("parse should succeed");
501 let atoms: Vec<_> = mol.atoms().collect();
502 assert_eq!(atoms[0].1.element, Element::C);
503 assert_eq!(atoms[1].1.element, Element::C);
504 assert_eq!(atoms[2].1.element, Element::O);
505 }
506
507 #[test]
508 fn test_parse_bond_types() {
509 let mol_str = "\
511test
512 chematic
513
514 8 4 0 0 0 0 0 0 0 0 0 V2000
515 0.0000 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0
516 1.0000 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0
517 2.0000 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0
518 3.0000 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0
519 4.0000 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0
520 5.0000 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0
521 6.0000 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0
522 7.0000 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0
523 1 2 1 0
524 3 4 2 0
525 5 6 3 0
526 7 8 4 0
527M END
528";
529 let (mol, _) = parse_mol(mol_str).expect("parse should succeed");
530 let bonds: Vec<_> = mol.bonds().collect();
531 assert_eq!(bonds[0].1.order, BondOrder::Single);
532 assert_eq!(bonds[1].1.order, BondOrder::Double);
533 assert_eq!(bonds[2].1.order, BondOrder::Triple);
534 assert_eq!(bonds[3].1.order, BondOrder::Aromatic);
535 }
536
537 #[test]
538 fn test_parse_query_bond_types_preserved() {
539 let mol_str = "\
540query_bonds
541 chematic
542
543 8 4 0 0 0 0 0 0 0 0 0 V2000
544 0.0000 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0
545 1.0000 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0
546 2.0000 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0
547 3.0000 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0
548 4.0000 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0
549 5.0000 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0
550 6.0000 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0
551 7.0000 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0
552 1 2 5 0
553 3 4 6 0
554 5 6 7 0
555 7 8 8 0
556M END
557";
558 let (mol, meta) = parse_mol(mol_str).expect("parse query bonds");
559 let bonds: Vec<_> = mol.bonds().collect();
560 assert_eq!(bonds[0].1.order, BondOrder::QuerySingleOrDouble);
561 assert_eq!(bonds[1].1.order, BondOrder::QuerySingleOrAromatic);
562 assert_eq!(bonds[2].1.order, BondOrder::QueryDoubleOrAromatic);
563 assert_eq!(bonds[3].1.order, BondOrder::QueryAny);
564
565 let written = write_mol(&mol, &meta);
566 assert!(written.contains(" 1 2 5 0"), "{written}");
567 assert!(written.contains(" 7 8 8 0"), "{written}");
568 }
569
570 #[test]
571 fn test_parse_charge() {
572 let mol_str = "\
574charged
575 chematic
576
577 1 0 0 0 0 0 0 0 0 0 0 V2000
578 0.0000 0.0000 0.0000 N 0 3 0 0 0 0 0 0 0 0 0 0
579M END
580";
581 let (mol, _) = parse_mol(mol_str).expect("parse should succeed");
582 assert_eq!(mol.atom(AtomIdx(0)).charge, 1);
583 }
584
585 #[test]
586 fn test_parse_negative_charge() {
587 let mol_str = "\
589negcharge
590 chematic
591
592 1 0 0 0 0 0 0 0 0 0 0 V2000
593 0.0000 0.0000 0.0000 O 0 5 0 0 0 0 0 0 0 0 0 0
594M END
595";
596 let (mol, _) = parse_mol(mol_str).expect("parse should succeed");
597 assert_eq!(mol.atom(AtomIdx(0)).charge, -1);
598 }
599
600 #[test]
601 fn test_round_trip() {
602 let (mol1, meta1) = parse_mol(ETHANOL_MOL).expect("first parse");
604 let written = write_mol(&mol1, &meta1);
605 let (mol2, _meta2) = parse_mol(&written).expect("second parse");
606 assert_eq!(mol1.atom_count(), mol2.atom_count());
607 assert_eq!(mol1.bond_count(), mol2.bond_count());
608 }
609
610 #[test]
611 fn test_round_trip_elements_preserved() {
612 let (mol1, meta1) = parse_mol(ETHANOL_MOL).expect("first parse");
613 let written = write_mol(&mol1, &meta1);
614 let (mol2, _) = parse_mol(&written).expect("second parse");
615 for ((_, a1), (_, a2)) in mol1.atoms().zip(mol2.atoms()) {
616 assert_eq!(a1.element, a2.element);
617 }
618 }
619
620 #[test]
621 fn test_error_missing_v2000() {
622 let bad = "\
623bad
624 prog
625
626 3 2 0 0 0 0 0 0 0 0 0 V3000
627 0.0000 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0
628M END
629";
630 assert!(matches!(
631 parse_mol(bad),
632 Err(MolParseError::InvalidCountLine { .. })
633 ));
634 }
635
636 #[test]
637 fn test_error_truncated_input() {
638 let bad = "\
640trunc
641 prog
642
643 3 0 0 0 0 0 0 0 0 0 0 V2000
644 0.0000 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0
645";
646 assert!(matches!(parse_mol(bad), Err(MolParseError::UnexpectedEnd)));
647 }
648
649 #[test]
650 fn test_error_invalid_counts_line() {
651 let bad = "\
653mol
654 prog
655
656 X Y
657M END
658";
659 assert!(matches!(
660 parse_mol(bad),
661 Err(MolParseError::InvalidCountLine { .. })
662 ));
663 }
664
665 #[test]
666 fn test_write_contains_m_end() {
667 let (mol, meta) = parse_mol(ETHANOL_MOL).expect("parse");
668 let written = write_mol(&mol, &meta);
669 assert!(written.contains("M END"));
670 }
671
672 #[test]
673 fn test_write_contains_v2000() {
674 let (mol, meta) = parse_mol(ETHANOL_MOL).expect("parse");
675 let written = write_mol(&mol, &meta);
676 assert!(written.contains("V2000"));
677 }
678
679 #[test]
680 fn test_parse_stereo_up_bond() {
681 let mol_str = "\n\n\n 2 1 0 0 0 0 999 V2000\n 0.0000 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0\n 1.0000 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0\n 1 2 1 1 0 0 0\nM END\n";
683 let (mol, _) = crate::parse_mol(mol_str).unwrap();
684 let bond = mol.bond(chematic_core::BondIdx(0));
685 assert_eq!(bond.order, chematic_core::BondOrder::Up);
686 }
687
688 #[test]
689 fn test_parse_stereo_down_bond() {
690 let mol_str = "\n\n\n 2 1 0 0 0 0 999 V2000\n 0.0000 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0\n 1.0000 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0\n 1 2 1 6 0 0 0\nM END\n";
691 let (mol, _) = crate::parse_mol(mol_str).unwrap();
692 let bond = mol.bond(chematic_core::BondIdx(0));
693 assert_eq!(bond.order, chematic_core::BondOrder::Down);
694 }
695
696 #[test]
697 fn test_molmetadata_builder() {
698 let meta = MolMetadata::default()
699 .with_name("aspirin")
700 .with_comment("test molecule");
701 assert_eq!(meta.name, "aspirin");
702 assert_eq!(meta.comment, "test molecule");
703 }
704
705 #[test]
706 fn test_molmetadata_with_name_roundtrip() {
707 use chematic_core::{Atom, BondOrder, Element, MoleculeBuilder};
709 let mut b = MoleculeBuilder::new();
710 let c1 = b.add_atom(Atom::new(Element::C));
711 let c2 = b.add_atom(Atom::new(Element::C));
712 b.add_bond(c1, c2, BondOrder::Single).unwrap();
713 let mol = b.build();
714
715 let meta = MolMetadata::default().with_name("acetic acid");
716 let molblock = crate::write_mol(&mol, &meta);
717 assert!(
718 molblock.starts_with("acetic acid"),
719 "MOL block must start with the molecule name"
720 );
721 }
722
723 #[test]
724 fn test_declared_max_atom_count_truncated_input_errors() {
725 let bad = "\
726max_atoms
727 chematic
728
729999 0 0 0 0 0 0 0 0 0 0 V2000
730";
731 assert!(matches!(parse_mol(bad), Err(MolParseError::UnexpectedEnd)));
732 }
733
734 #[test]
735 fn test_declared_large_bond_count_truncated_input_errors() {
736 let bad = "\
737many_bonds
738 chematic
739
740 1 999 0 0 0 0 0 0 0 0 0 V2000
741 0.0000 0.0000 0.0000 C 0 0 0 0 0 0 0 0 0 0 0 0
742";
743 assert!(matches!(parse_mol(bad), Err(MolParseError::UnexpectedEnd)));
744 }
745}