com.actelion.research.chem.phesa.MolecularVolume Maven / Gradle / Ivy
Go to download
Show more of this group Show more artifacts with this name
Show all versions of openchemlib Show documentation
Show all versions of openchemlib Show documentation
Open Source Chemistry Library
package com.actelion.research.chem.phesa;
import com.actelion.research.calc.Matrix;
import com.actelion.research.chem.Coordinates;
import com.actelion.research.chem.StereoMolecule;
import com.actelion.research.chem.alignment3d.transformation.ExponentialMap;
import com.actelion.research.chem.alignment3d.transformation.Quaternion;
import com.actelion.research.chem.alignment3d.transformation.Transformation;
import com.actelion.research.chem.conf.Conformer;
import com.actelion.research.chem.phesa.pharmacophore.IonizableGroupDetector;
import com.actelion.research.chem.phesa.pharmacophore.PharmacophoreCalculator;
import com.actelion.research.chem.phesa.pharmacophore.pp.ExitVectorPoint;
import com.actelion.research.chem.phesa.pharmacophore.pp.IPharmacophorePoint;
import com.actelion.research.chem.phesa.pharmacophore.pp.PPGaussian;
import java.util.ArrayList;
import java.util.List;
import java.util.stream.Collectors;
import com.actelion.research.util.EncoderFloatingPointNumbers;
/**
* @version: 1.0, February 2018
* Author: J. Wahl
* class to approximate the volume of a molecule as a sum of atom-centered Gaussians, as introduced by Grant and Pickup, J. Phys. Chem. 1995, 99, 3503-3510
* no higher order terms (atom-atom overlaps) are calculated to reduce computational costs
*/
public class MolecularVolume extends ShapeVolume{
static public final double p = 2.82842712475; // height of Atomic Gaussian, 2*sqrt(2), commonly used in the literature: Haque and Pande, DOI 10.1002/jcc.11307
static public final double alpha_pref = 2.41798793102; // taken from DOI 10.1002/jcc.11307
private ArrayList volumeGaussians; //exclusion and inclusion spheres
private ArrayList hydrogens;
public MolecularVolume(List atomicGaussiansInp,List ppGaussiansInp,List volGaussians,
List hydrogenCoords) {
this.atomicGaussians = new ArrayList();
for(AtomicGaussian ag : atomicGaussiansInp) {
atomicGaussians.add(new AtomicGaussian(ag));
}
this.ppGaussians = new ArrayList();
for(PPGaussian pg : ppGaussiansInp) {
ppGaussians.add(new PPGaussian(pg));
}
this.hydrogens = new ArrayList();
for(Coordinates hydrogen : hydrogenCoords) {
this.hydrogens.add(hydrogen);
}
this.volumeGaussians = new ArrayList();
for(VolumeGaussian eg : volGaussians) {
this.volumeGaussians.add(new VolumeGaussian(eg));
}
this.calcCOM();
}
public void updateCOM() {
this.calcCOM();
}
private void updateAtomIndeces(List gaussians,int[] map) {
for(Gaussian3D gaussian:gaussians)
gaussian.updateAtomIndeces(map);
}
public void updateAtomIndeces(int[] map) {
updateAtomIndeces(ppGaussians,map);
updateAtomIndeces(atomicGaussians,map);
updateAtomIndeces(volumeGaussians,map);
}
public MolecularVolume(StereoMolecule mol) {
this.hydrogens = new ArrayList();
this.volumeGaussians = new ArrayList();
this.calc(mol);
this.calcPPVolume(mol);
this.calcCOM();
}
public MolecularVolume(MolecularVolume original, Conformer conf) {
this(original);
update(conf);
}
/**
* calculate the Overlap of the two molecular volumes as a function a transform vector that is applied to the query molecule
* overlap Volume of two molecular Volumes formulated as a summed overlap of atomic Gaussians
* taken from Grant, Gallardo, Pickup, Journal of Computational Chemistry, 17, 1653-1666, 1996
* returns a double[2]: the first double is the total overlap, whereas the second value is the specific
* contribution of additional volume gaussians (inclusion, exclusion)
* @param transform
* @return
*/
public MolecularVolume(MolecularVolume original) {
this.atomicGaussians = new ArrayList();
this.ppGaussians = new ArrayList();
this.volumeGaussians = new ArrayList();
for(AtomicGaussian ag : original.getAtomicGaussians()) {
this.atomicGaussians.add(new AtomicGaussian(ag));
}
for(PPGaussian pg : original.getPPGaussians()) {
this.ppGaussians.add(new PPGaussian(pg));
}
this.hydrogens = new ArrayList();
for(Coordinates hydrogen : original.hydrogens) {
this.hydrogens.add(new Coordinates(hydrogen));
}
for(VolumeGaussian eg : original.volumeGaussians) {
this.volumeGaussians.add(new VolumeGaussian(eg));
}
this.com = new Coordinates(original.com);
}
/**
* calculate the self-overlap of the base molecule
* @return
*/
@Override
public double getSelfAtomOverlap(){
double Vtot = 0.0;
for(AtomicGaussian at:atomicGaussians){
for(AtomicGaussian at2:atomicGaussians){
Vtot += at.getVolumeOverlap(at2);
}
for(VolumeGaussian vg : volumeGaussians) {
if(vg.getRole()!=VolumeGaussian.INCLUSION)
continue;
Vtot += vg.getRole()*at.getVolumeOverlap(vg);
}
}
for(VolumeGaussian vg:volumeGaussians){
if(vg.getRole()!=VolumeGaussian.INCLUSION)
continue;
for(VolumeGaussian vg2 : volumeGaussians) {
//only consider self-overlap of inclusion spheres
if(vg2.getRole()!=VolumeGaussian.INCLUSION)
continue;
Vtot += vg2.getVolumeOverlap(vg);
}
}
return Vtot;
}
public double[] getTotalAtomOverlap(double[] transform, MolecularVolume fitVol){
double[] result = new double[2];
ExponentialMap eMap = new ExponentialMap(transform[0],transform[1],transform[2]);
double Vtot = 0.0;
double Vvol = 0.0;
Coordinates com = fitVol.getCOM();
double[][] rotMatrix = eMap.toQuaternion().getRotMatrix().getArray();
List fitGaussians = fitVol.atomicGaussians;
Coordinates[] fitCenterModCoords = new Coordinates[fitGaussians.size()];
for(int k=0;k();
int nrOfAtoms = mol.getAllAtoms();
for (int i=0;i();
List ppPoints = new ArrayList();
IonizableGroupDetector detector = new IonizableGroupDetector(mol);
ppPoints.addAll(detector.detect());
ppPoints.addAll(PharmacophoreCalculator.getPharmacophorePoints(mol));
for(IPharmacophorePoint ppPoint : ppPoints )
ppGaussians.add(new PPGaussian(6,ppPoint));
}
/**
* calculates volume weighted center of mass of the molecular Volume
*/
@Override
public void calcCOM(){
double volume = 0.0;
double comX = 0.0;
double comY = 0.0;
double comZ = 0.0;
for(AtomicGaussian atGauss : this.atomicGaussians){
volume += atGauss.getVolume();
comX += atGauss.getCenter().x*atGauss.getVolume();
comY += atGauss.getCenter().y*atGauss.getVolume();
comZ += atGauss.getCenter().z*atGauss.getVolume();
}
for(VolumeGaussian volGauss : this.volumeGaussians){
volume += volGauss.getRole()*volGauss.getVolume();
comX += volGauss.getRole()*volGauss.getCenter().x*volGauss.getVolume();
comY += volGauss.getRole()*volGauss.getCenter().y*volGauss.getVolume();
comZ += volGauss.getRole()*volGauss.getCenter().z*volGauss.getVolume();
}
comX = comX/volume;
comY = comY/volume;
comZ = comZ/volume;
this.com = new Coordinates(comX,comY,comZ);
}
@Override
public Matrix getCovarianceMatrix() {
Matrix massMatrix = new Matrix(3,3);
double volume = 0.0;
for (AtomicGaussian ag : atomicGaussians){
volume += ag.getVolume();
double value = ag.getVolume()*ag.getCenter().x*ag.getCenter().x;
massMatrix.addToElement(0,0,value);
value = ag.getVolume()*ag.getCenter().x*ag.getCenter().y;
massMatrix.addToElement(0,1,value);
value = ag.getVolume()*ag.getCenter().x*ag.getCenter().z;
massMatrix.addToElement(0,2,value);
value = ag.getVolume()*ag.getCenter().y*ag.getCenter().y;
massMatrix.addToElement(1,1,value);
value = ag.getVolume()*ag.getCenter().y*ag.getCenter().z;
massMatrix.addToElement(1,2,value);
value = ag.getVolume()*ag.getCenter().z*ag.getCenter().z;
massMatrix.addToElement(2,2,value);
}
for (VolumeGaussian vg : volumeGaussians){
volume += vg.getRole()*vg.getVolume();
double value = vg.getRole()*vg.getVolume()*vg.getCenter().x*vg.getCenter().x;
massMatrix.addToElement(0,0,value);
value = vg.getRole()*vg.getVolume()*vg.getCenter().x*vg.getCenter().y;
massMatrix.addToElement(0,1,value);
value = vg.getRole()*vg.getVolume()*vg.getCenter().x*vg.getCenter().z;
massMatrix.addToElement(0,2,value);
value = vg.getRole()*vg.getVolume()*vg.getCenter().y*vg.getCenter().y;
massMatrix.addToElement(1,1,value);
value = vg.getRole()*vg.getVolume()*vg.getCenter().y*vg.getCenter().z;
massMatrix.addToElement(1,2,value);
value = vg.getRole()*vg.getVolume()*vg.getCenter().z*vg.getCenter().z;
massMatrix.addToElement(2,2,value);
}
massMatrix.set(0,0,massMatrix.get(0,0)/volume);
massMatrix.set(0,1,massMatrix.get(0,1)/volume);
massMatrix.set(0,2,massMatrix.get(0,2)/volume);
massMatrix.set(1,1,massMatrix.get(1,1)/volume);
massMatrix.set(1,2,massMatrix.get(1,2)/volume);
massMatrix.set(2,2,massMatrix.get(2,2)/volume);
massMatrix.set(1,0,massMatrix.get(0,1));
massMatrix.set(2,0,massMatrix.get(0,2));
massMatrix.set(2,1,massMatrix.get(1,2));
return massMatrix;
}
public ArrayList getVolumeGaussians() {
return this.volumeGaussians;
}
@Override
protected void transformGaussians(Transformation transform) {
super.transformGaussians(transform);
transformGaussians(volumeGaussians,transform);
}
public ArrayList getHydrogens() {
return this.hydrogens;
}
private void updateHydrogens(StereoMolecule mol) {
int h = 0;
for(int i = mol.getAtoms();i referenceAtomicGaussians = reference.getAtomicGaussians();
List referencePPGaussians = reference.getPPGaussians();
List referenceVolGaussians = reference.getVolumeGaussians();
String[] splitString = string.split(" ");
double[] atomicGaussiansCoords = EncoderFloatingPointNumbers.decode(splitString[0]);
double[] ppGaussiansCoords = EncoderFloatingPointNumbers.decode(splitString[1]);
double[] ppGaussiansDirectionalities = EncoderFloatingPointNumbers.decode(splitString[2]);
double[] volumeGaussiansRefCoords = EncoderFloatingPointNumbers.decode(splitString[3]);
double[] volumeGaussiansShiftCoords = EncoderFloatingPointNumbers.decode(splitString[4]);
double[] hydrogensCoords = EncoderFloatingPointNumbers.decode(splitString[5]);
ArrayList atomicGaussians = new ArrayList();
ArrayList ppGaussians = new ArrayList();
ArrayList volumeGaussians = new ArrayList();
ArrayList hydrogens = new ArrayList();
int nrOfAtomicGaussians = atomicGaussiansCoords.length/3;
int nrOfHydrogens = hydrogensCoords.length/3;
int nrOfPPGaussians = ppGaussiansCoords.length/3;
int nrOfVolumeGaussians = volumeGaussiansRefCoords.length/3;
for(int i=0;i atomicGaussians = new ArrayList();
List ppGaussians = new ArrayList();
List volumeGaussians = new ArrayList();
List hydrogens = new ArrayList();
for(int i=firstIndex;i