Skip to content

Commit

Permalink
moving generation of atom tags to chem.soft.-agnostic code
Browse files Browse the repository at this point in the history
  • Loading branch information
marco-foscato committed Nov 21, 2023
1 parent e213602 commit dc6e85d
Show file tree
Hide file tree
Showing 13 changed files with 170 additions and 48 deletions.
30 changes: 9 additions & 21 deletions src/main/java/autocompchem/atom/AtomUtils.java
Original file line number Diff line number Diff line change
Expand Up @@ -345,15 +345,7 @@ public static String getSymbolOrLabel(IAtom atm)
}
else if (atm instanceof PseudoAtom)
{
if (isPsaudoAtmWithLabel(atm,AtomConstants.DUMMYATMLABEL))
{
res = AtomConstants.DUMMYATMLABEL;
}
else if (isPsaudoAtmWithLabel(atm,
AtomConstants.ATTACHMENTPOINTLABEL))
{
res = AtomConstants.ATTACHMENTPOINTLABEL;
}
res = ((PseudoAtom) atm).getLabel();
}
else if (atm instanceof Atom)
{
Expand Down Expand Up @@ -418,18 +410,14 @@ public static Integer countExplicitHydrogens(IAtom atm, IAtomContainer mol)
public static Point3d getCoords3d(IAtom atom)
{
Point3d p3d = new Point3d();
try {
Point2d atp2d = new Point2d();
atp2d = atom.getPoint2d();
p3d.x = atp2d.x;
p3d.y = atp2d.y;
p3d.z = 0.0000;
} catch (Throwable t) {
Point3d atp3d = new Point3d();
atp3d = atom.getPoint3d();
p3d.x = atp3d.x;
p3d.y = atp3d.y;
p3d.z = atp3d.z;
if (atom.getPoint3d()!=null)
{
Point3d atp3d = atom.getPoint3d();
p3d = new Point3d(atp3d.x, atp3d.y, atp3d.z);
} else if (atom.getPoint2d()!=null)
{
Point2d atp2d = atom.getPoint2d();
p3d = new Point3d(atp2d.x, atp2d.y, 0);
}
return p3d;
}
Expand Down
23 changes: 20 additions & 3 deletions src/main/java/autocompchem/chemsoftware/ChemSoftInputWriter.java
Original file line number Diff line number Diff line change
Expand Up @@ -515,8 +515,12 @@ private void produceSingleJobInputFiles(List<IAtomContainer> mols,
molSpecJob.setParameter(ChemSoftConstants.PAROUTFILEROOT,
outFileNameRoot, true);

// Add atom coordinates to the so-far molecule-agnostic job
setChemicalSystem(molSpecJob, mols);
// Add atom coordinates to the so-far possible molecule-agnostic job
if (useAtomTags)
setChemicalSystem(molSpecJob, makeAtomContainersWithAtomTags(mols));
else
setChemicalSystem(molSpecJob, mols);

// We keep a copy of the agnostic chemical system definition in the job
// parameters
setChemicalSystemAsJobParam(molSpecJob, mols);
Expand All @@ -543,7 +547,7 @@ private void produceSingleJobInputFiles(List<IAtomContainer> mols,
}
}

// These calls take care also of the sub-jobs/directives
// This call takes care also of the sub-jobs/directives
molSpecJob.processDirectives(mols, this.getMyJob());

// Ensure a value of charge and spin has been defined
Expand All @@ -570,6 +574,19 @@ private void produceSingleJobInputFiles(List<IAtomContainer> mols,
}
}

//------------------------------------------------------------------------------

private List<IAtomContainer> makeAtomContainersWithAtomTags(
List<IAtomContainer> mols)
{
List<IAtomContainer> molsWithAtomTags = new ArrayList<IAtomContainer>();
for (IAtomContainer iac : mols)
{
molsWithAtomTags.add(MolecularUtils.makeSimpleCopyWithAtomTags(iac));
}
return molsWithAtomTags;
}

//------------------------------------------------------------------------------

private void setChemicalSystemAsJobParam(CompChemJob job,
Expand Down
8 changes: 8 additions & 0 deletions src/main/java/autocompchem/chemsoftware/Directive.java
Original file line number Diff line number Diff line change
Expand Up @@ -47,6 +47,7 @@
import autocompchem.modeling.basisset.BasisSetGenerator;
import autocompchem.modeling.constraints.ConstraintsGenerator;
import autocompchem.modeling.constraints.ConstraintsSet;
import autocompchem.molecule.MolecularUtils;
import autocompchem.molecule.conformation.ConformationalSpace;
import autocompchem.molecule.conformation.ConformationalSpaceGenerator;
import autocompchem.molecule.intcoords.zmatrix.ZMatrix;
Expand Down Expand Up @@ -1231,6 +1232,8 @@ private void performACCTask(List<IAtomContainer> mols, ParameterStorage params,
ChemSoftConstants.PARMULTIGEOMID).getValueAsString());
}

boolean useAtomTags = params.contains(ChemSoftConstants.PARUSEATMTAGS);

switch (coordsType)
{
case ZMAT:
Expand All @@ -1244,6 +1247,7 @@ private void performACCTask(List<IAtomContainer> mols, ParameterStorage params,
Worker w = WorkerFactory.createWorker(zmatMakerTask,
masterJob);
ZMatrixHandler zmh = (ZMatrixHandler) w;
//TODO-gg make it use atom tags upon request
ZMatrix zmat = zmh.makeZMatrix();
if (params.contains(ZMatrixConstants.SELECTORMODE))
{
Expand Down Expand Up @@ -1290,6 +1294,10 @@ private void performACCTask(List<IAtomContainer> mols, ParameterStorage params,
default:
{
IAtomContainer mol = mols.get(geometryId);
if (useAtomTags)
{
mol = MolecularUtils.makeSimpleCopyWithAtomTags(mol);
}
((IValueContainer) dirComp).setValue(mol);
break;
}
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -272,22 +272,10 @@ private ArrayList<String> getTextForInput(Directive d)
break;

case IATOMCONTAINER:
boolean useTags = false;
if (data.getTaskParams()!=null
&& data.getTaskParams().contains(
ChemSoftConstants.PARUSEATMTAGS))
{
useTags = true;
}
IAtomContainer mol = (IAtomContainer) data.getValue();
for (IAtom atm : mol.atoms())
{
String atmId = AtomUtils.getSymbolOrLabel(atm);
if (useTags)
{
// Convention is to use 1-based indexing here
atmId = atmId + (mol.indexOf(atm)+1);
}
Point3d p3d = AtomUtils.getCoords3d(atm);
ddLines.add(String.format(Locale.ENGLISH," %3s", atmId)
+ String.format(Locale.ENGLISH," %10.6f",p3d.x)
Expand Down
54 changes: 53 additions & 1 deletion src/main/java/autocompchem/molecule/MolecularUtils.java
Original file line number Diff line number Diff line change
Expand Up @@ -30,13 +30,19 @@
import javax.vecmath.Point3d;
import javax.vecmath.Vector3d;

import org.openscience.cdk.Atom;
import org.openscience.cdk.Bond;
import org.openscience.cdk.CDKConstants;
import org.openscience.cdk.DefaultChemObjectBuilder;
import org.openscience.cdk.PseudoAtom;
import org.openscience.cdk.aromaticity.Kekulization;
import org.openscience.cdk.exception.CDKException;
import org.openscience.cdk.interfaces.IAtom;
import org.openscience.cdk.interfaces.IAtomContainer;
import org.openscience.cdk.interfaces.IBond;
import org.openscience.cdk.interfaces.IChemObjectBuilder;
import org.openscience.cdk.interfaces.IBond.Order;
import org.openscience.cdk.tools.periodictable.PeriodicTable;

import autocompchem.atom.AtomUtils;
import autocompchem.geometry.DistanceMatrix;
Expand All @@ -58,6 +64,9 @@ public class MolecularUtils
//TODO move to constants
@SuppressWarnings("unused")
private static double linearBendThld = 175.0;

private static IChemObjectBuilder chemBuilder =
DefaultChemObjectBuilder.getInstance();


//------------------------------------------------------------------------------
Expand Down Expand Up @@ -807,7 +816,50 @@ public static void ensureNoUnsetBondOrders(IAtomContainer mol)
}
}
}


//------------------------------------------------------------------------------

/**
* Makes a simplified copy of the given {@link IAtomContainer} considering
* only the list of atoms as elements/labels and points in 3D, set of bonds,
* and container title/name,
* without any further property of atoms, bonds, or of the container, and
* sets the atom label to be an atom tag, i.e., a string that allows to
* identify the atom (e.g., "C3"). The tag is generated by concatenating
* the elemental symbol (or initial {@link PseudoAtom} label) with the
* 1-based index of the atom in the given container.
* @param iac the original atom container to copy and decorate with atom
* tags.
* @return the copy of the original container with the atom tags.
*/
public static IAtomContainer makeSimpleCopyWithAtomTags(IAtomContainer iac)
{
IAtomContainer iacWithTags = chemBuilder.newAtomContainer();
int iAtm = 0;
for (IAtom origAtm : iac.atoms())
{
iAtm++;
String atomTag = AtomUtils.getSymbolOrLabel(origAtm) + iAtm;
IAtom atmWithTag = new PseudoAtom(atomTag);
Point3d p3d = AtomUtils.getCoords3d(origAtm);
atmWithTag.setPoint3d(new Point3d(p3d.x, p3d.y, p3d.z));
iacWithTags.addAtom(atmWithTag);
}
for (IBond origBnd : iac.bonds())
{
IAtom[] atomsInNewBond = new IAtom[origBnd.getAtomCount()];
for (int i=0; i<origBnd.getAtomCount(); i++)
{
atomsInNewBond[i] = iacWithTags.getAtom(
iac.indexOf(origBnd.getAtom(i)));
}
IBond newBond = new Bond(atomsInNewBond, origBnd.getOrder());
iacWithTags.addBond(newBond);
}
iacWithTags.setTitle(MolecularUtils.getNameIDOrNull(iac));
return iacWithTags;
}

//------------------------------------------------------------------------------

}
37 changes: 37 additions & 0 deletions src/test/java/autocompchem/atom/AtomUtilsTest.java
Original file line number Diff line number Diff line change
Expand Up @@ -24,10 +24,17 @@

import static org.junit.jupiter.api.Assertions.assertTrue;

import javax.vecmath.Point2d;
import javax.vecmath.Point3d;

import org.junit.jupiter.api.Test;
import org.openscience.cdk.Atom;
import org.openscience.cdk.DefaultChemObjectBuilder;
import org.openscience.cdk.PseudoAtom;
import org.openscience.cdk.interfaces.IAtom;
import org.openscience.cdk.interfaces.IChemObjectBuilder;

import autocompchem.molecule.MolecularUtils;

/**
* Unit Test for atom utilities
Expand All @@ -38,6 +45,9 @@
public class AtomUtilsTest
{

private static IChemObjectBuilder chemBuilder =
DefaultChemObjectBuilder.getInstance();

//------------------------------------------------------------------------------

@Test
Expand Down Expand Up @@ -81,6 +91,33 @@ public void testGetSymbolOrLabel() throws Exception
"Return label for attachment points");
}

//------------------------------------------------------------------------------

@Test
public void testGetCoords3d() throws Exception
{
IAtom atm = chemBuilder.newAtom();
assertTrue(closeEnough(AtomUtils.getCoords3d(atm), new Point3d(0,0,0)));

atm.setPoint2d(new Point2d(1.2, 2.3));
assertTrue(closeEnough(AtomUtils.getCoords3d(atm),
new Point3d(1.2, 2.3,0)));

atm.setPoint3d(new Point3d(1.2, 2.3, -3.4));
assertTrue(closeEnough(AtomUtils.getCoords3d(atm),
new Point3d(1.2, 2.3, -3.4)));
}

//------------------------------------------------------------------------------

private boolean closeEnough(Point3d pA, Point3d pB)
{
if (pA==null | pB==null)
return false;

return pA.distance(pB) < 0.0002;
}

//------------------------------------------------------------------------------

}
43 changes: 38 additions & 5 deletions src/test/java/autocompchem/molecule/MolecularUtilsTest.java
Original file line number Diff line number Diff line change
@@ -1,5 +1,6 @@
package autocompchem.molecule;

import static org.junit.jupiter.api.Assertions.assertEquals;
import static org.junit.jupiter.api.Assertions.assertFalse;

/*
Expand All @@ -25,10 +26,13 @@

import org.junit.jupiter.api.Test;
import org.openscience.cdk.Atom;
import org.openscience.cdk.AtomContainer;
import org.openscience.cdk.DefaultChemObjectBuilder;
import org.openscience.cdk.PseudoAtom;
import org.openscience.cdk.interfaces.IAtom;
import org.openscience.cdk.interfaces.IAtomContainer;
import org.openscience.cdk.interfaces.IChemObjectBuilder;

import autocompchem.atom.AtomUtils;
import autocompchem.geometry.DistanceMatrix;

/**
Expand All @@ -39,13 +43,43 @@

public class MolecularUtilsTest
{

private static IChemObjectBuilder chemBuilder =
DefaultChemObjectBuilder.getInstance();

//------------------------------------------------------------------------------

@Test
public void testMakeSimpleCopyWithAtomTags() throws Exception
{
IAtomContainer mol = chemBuilder.newAtomContainer();
mol.addAtom(new Atom("H"));
mol.addAtom(new Atom("C"));
mol.addAtom(new Atom("O"));
mol.addAtom(new Atom("Ru"));
mol.addAtom(new PseudoAtom("Du"));
mol.addAtom(new PseudoAtom("Xx"));

IAtomContainer tagged = MolecularUtils.makeSimpleCopyWithAtomTags(mol);

assertEquals(mol.getAtomCount(),tagged.getAtomCount());
for (int iAtm=0; iAtm<mol.getAtomCount(); iAtm++)
{
IAtom atm = tagged.getAtom(iAtm);
IAtom origAtm = mol.getAtom(iAtm);
// NB: so far there is only one format for atom tags, so we can
// hard-code the expected tag.
String expectedTag = AtomUtils.getSymbolOrLabel(origAtm)+(iAtm+1);
assertEquals(expectedTag, AtomUtils.getSymbolOrLabel(atm));
}
}

//------------------------------------------------------------------------------

@Test
public void testElementalSymbols() throws Exception
{
IAtomContainer mol = new AtomContainer();
IAtomContainer mol = chemBuilder.newAtomContainer();
mol.addAtom(new Atom("H"));
mol.addAtom(new Atom("C"));
mol.addAtom(new Atom("O"));
Expand All @@ -72,7 +106,7 @@ public void testElementalSymbols() throws Exception
@Test
public void testMinInterelementalBondDistance() throws Exception
{
IAtomContainer mol = new AtomContainer();
IAtomContainer mol = chemBuilder.newAtomContainer();
mol.addAtom(new Atom("C",new Point3d(0,0,0.0)));
mol.addAtom(new Atom("C",new Point3d(1.0,0,0)));
mol.addAtom(new Atom("C",new Point3d(2.0,0,0)));
Expand All @@ -92,8 +126,7 @@ public void testMinInterelementalBondDistance() throws Exception
@Test
public void testInteratormicDistanceMatrix() throws Exception
{

IAtomContainer mol = new AtomContainer();
IAtomContainer mol = chemBuilder.newAtomContainer();
mol.addAtom(new Atom("C",new Point3d(0,0,0.0)));
mol.addAtom(new Atom("C",new Point3d(0,3.0,0)));
mol.addAtom(new Atom("C",new Point3d(4.0,0,0)));
Expand Down
2 changes: 1 addition & 1 deletion test/t87.check
Original file line number Diff line number Diff line change
Expand Up @@ -16,7 +16,7 @@ if [ $n -ne 2 ] ; then echo NOT Passes ERROR 3 ; exit -1 ; fi
n=0;n=$(grep -c -i "Ru.*nelec 28" t87-mol.nw)
if [ $n -ne 1 ] ; then echo NOT Passes ERROR 4 ; exit -1 ; fi

n=0;n=$(grep -c -i "Ru.*\\s*[0-9\-]*\.[0-9]*\\s*[0-9\-]*\.[0-9]*\\s*[0-9\-]*\.[0-9]*" t87-mol.nw)
n=0;n=$(grep -c -i "Ru1.*\\s*2.7141[0-9]*\\s*-1\.3069[0-9]*\\s*-0\.6133[0-9]*" t87-mol.nw)
if [ $n -ne 1 ] ; then echo NOT Passes ERROR 5 ; exit -1 ; fi

n=0;n=$(grep -c -i "fix atom " t87-mol.nw)
Expand Down
1 change: 1 addition & 0 deletions test/t87.params
Original file line number Diff line number Diff line change
Expand Up @@ -4,4 +4,5 @@
VERBOSITY: 2
TASK: PrepareInputNWChem
INPUTGEOMETRIESFILE: ../t87-mol.sdf
USEATOMTAGS:
JOBDETAILSFILE: ../t87-jobDetails.json
Loading

0 comments on commit dc6e85d

Please sign in to comment.