use of org.openscience.cdk.charges.Polarizability in project cdk by cdk.
the class BCUTDescriptor method calculate.
/**
* Calculates the three classes of BCUT descriptors.
*
* @param container Parameter is the atom container.
* @return An ArrayList containing the descriptors. The default is to return
* all calculated eigenvalues of the Burden matrices in the order described
* above. If a parameter list was supplied, then only the specified number
* of highest and lowest eigenvalues (for each class of BCUT) will be returned.
*/
@Override
public DescriptorValue calculate(IAtomContainer container) {
int counter;
IAtomContainer molecule;
try {
molecule = container.clone();
} catch (CloneNotSupportedException e) {
logger.debug("Error during clone");
return getDummyDescriptorValue(new CDKException("Error occurred during clone " + e));
}
// add H's in case they're not present
try {
AtomContainerManipulator.percieveAtomTypesAndConfigureAtoms(molecule);
CDKHydrogenAdder hAdder = CDKHydrogenAdder.getInstance(molecule.getBuilder());
hAdder.addImplicitHydrogens(molecule);
AtomContainerManipulator.convertImplicitToExplicitHydrogens(molecule);
} catch (Exception e) {
return getDummyDescriptorValue(new CDKException("Could not add hydrogens: " + e.getMessage(), e));
}
// do aromaticity detecttion for calculating polarizability later on
if (this.checkAromaticity) {
try {
AtomContainerManipulator.percieveAtomTypesAndConfigureAtoms(molecule);
} catch (CDKException e) {
return getDummyDescriptorValue(new CDKException("Error in atom typing: " + e.getMessage(), e));
}
try {
Aromaticity.cdkLegacy().apply(molecule);
} catch (CDKException e) {
return getDummyDescriptorValue(new CDKException("Error in aromaticity perception: " + e.getMessage()));
}
}
// find number of heavy atoms
int nheavy = 0;
for (int i = 0; i < molecule.getAtomCount(); i++) {
if (molecule.getAtom(i).getAtomicNumber() != IElement.H)
nheavy++;
}
if (nheavy == 0)
return getDummyDescriptorValue(new CDKException("No heavy atoms in the molecule"));
double[] diagvalue = new double[nheavy];
// get atomic mass weighted BCUT
counter = 0;
try {
for (int i = 0; i < molecule.getAtomCount(); i++) {
if (molecule.getAtom(i).getAtomicNumber() == IElement.H)
continue;
diagvalue[counter] = Isotopes.getInstance().getMajorIsotope(molecule.getAtom(i).getSymbol()).getExactMass();
counter++;
}
} catch (Exception e) {
return getDummyDescriptorValue(new CDKException("Could not calculate weight: " + e.getMessage(), e));
}
double[][] burdenMatrix = BurdenMatrix.evalMatrix(molecule, diagvalue);
if (hasUndefined(burdenMatrix))
return getDummyDescriptorValue(new CDKException("Burden matrix has undefined values"));
Matrix matrix = new Matrix(burdenMatrix);
EigenvalueDecomposition eigenDecomposition = new EigenvalueDecomposition(matrix);
double[] eval1 = eigenDecomposition.getRealEigenvalues();
// get charge weighted BCUT
LonePairElectronChecker lpcheck = new LonePairElectronChecker();
GasteigerMarsiliPartialCharges peoe;
try {
lpcheck.saturate(molecule);
double[] charges = new double[molecule.getAtomCount()];
// pepe = new GasteigerPEPEPartialCharges();
// pepe.calculateCharges(molecule);
// for (int i = 0; i < molecule.getAtomCount(); i++) charges[i] = molecule.getAtom(i).getCharge();
peoe = new GasteigerMarsiliPartialCharges();
peoe.assignGasteigerMarsiliSigmaPartialCharges(molecule, true);
for (int i = 0; i < molecule.getAtomCount(); i++) charges[i] += molecule.getAtom(i).getCharge();
for (int i = 0; i < molecule.getAtomCount(); i++) {
molecule.getAtom(i).setCharge(charges[i]);
}
} catch (Exception e) {
return getDummyDescriptorValue(new CDKException("Could not calculate partial charges: " + e.getMessage(), e));
}
counter = 0;
for (int i = 0; i < molecule.getAtomCount(); i++) {
if (molecule.getAtom(i).getAtomicNumber() == IElement.H)
continue;
diagvalue[counter] = molecule.getAtom(i).getCharge();
counter++;
}
burdenMatrix = BurdenMatrix.evalMatrix(molecule, diagvalue);
if (hasUndefined(burdenMatrix))
return getDummyDescriptorValue(new CDKException("Burden matrix has undefined values"));
matrix = new Matrix(burdenMatrix);
eigenDecomposition = new EigenvalueDecomposition(matrix);
double[] eval2 = eigenDecomposition.getRealEigenvalues();
int[][] topoDistance = PathTools.computeFloydAPSP(AdjacencyMatrix.getMatrix(molecule));
// get polarizability weighted BCUT
Polarizability pol = new Polarizability();
counter = 0;
for (int i = 0; i < molecule.getAtomCount(); i++) {
if (molecule.getAtom(i).getAtomicNumber() == IElement.H)
continue;
diagvalue[counter] = pol.calculateGHEffectiveAtomPolarizability(molecule, molecule.getAtom(i), false, topoDistance);
counter++;
}
burdenMatrix = BurdenMatrix.evalMatrix(molecule, diagvalue);
if (hasUndefined(burdenMatrix))
return getDummyDescriptorValue(new CDKException("Burden matrix has undefined values"));
matrix = new Matrix(burdenMatrix);
eigenDecomposition = new EigenvalueDecomposition(matrix);
double[] eval3 = eigenDecomposition.getRealEigenvalues();
// return only the n highest & lowest eigenvalues
int lnlow, lnhigh, enlow, enhigh;
if (nlow > nheavy) {
lnlow = nheavy;
enlow = nlow - nheavy;
} else {
lnlow = nlow;
enlow = 0;
}
if (nhigh > nheavy) {
lnhigh = nheavy;
enhigh = nhigh - nheavy;
} else {
lnhigh = nhigh;
enhigh = 0;
}
DoubleArrayResult retval = new DoubleArrayResult((lnlow + enlow + lnhigh + enhigh) * 3);
for (int i = 0; i < lnlow; i++) retval.add(eval1[i]);
for (int i = 0; i < enlow; i++) retval.add(Double.NaN);
for (int i = 0; i < lnhigh; i++) retval.add(eval1[eval1.length - i - 1]);
for (int i = 0; i < enhigh; i++) retval.add(Double.NaN);
for (int i = 0; i < lnlow; i++) retval.add(eval2[i]);
for (int i = 0; i < enlow; i++) retval.add(Double.NaN);
for (int i = 0; i < lnhigh; i++) retval.add(eval2[eval2.length - i - 1]);
for (int i = 0; i < enhigh; i++) retval.add(Double.NaN);
for (int i = 0; i < lnlow; i++) retval.add(eval3[i]);
for (int i = 0; i < enlow; i++) retval.add(Double.NaN);
for (int i = 0; i < lnhigh; i++) retval.add(eval3[eval3.length - i - 1]);
for (int i = 0; i < enhigh; i++) retval.add(Double.NaN);
return new DescriptorValue(getSpecification(), getParameterNames(), getParameters(), retval, getDescriptorNames());
}
use of org.openscience.cdk.charges.Polarizability in project cdk by cdk.
the class IonizationPotentialTool method getQSARs.
/**
* Get the results of 7 qsar descriptors been applied. They are:
* Electronegativity,
* GasteigerMarsiliPartialCharges,
* GasteigerPEPEPartialCharges,
* Polarizability,
* StabilizationCharge,
* Number of Atom in resonance
* if the container in resonance is aromatic.
*
* @param container The IAtomContainer which contain the IAtom
* @param atom The IAtom to calculate
* @return An Array containing the results
*/
public static double[] getQSARs(IAtomContainer container, IAtom atom) throws CDKException {
Electronegativity electronegativity = new Electronegativity();
PiElectronegativity pielectronegativity = new PiElectronegativity();
GasteigerMarsiliPartialCharges peoe = new GasteigerMarsiliPartialCharges();
GasteigerPEPEPartialCharges pepe = new GasteigerPEPEPartialCharges();
Polarizability pol = new Polarizability();
StabilizationCharges stabil = new StabilizationCharges();
StructureResonanceGenerator gRI = new StructureResonanceGenerator();
IAtomContainer product = initiateIonization(container, atom);
double[] results = new double[8];
// sigmaElectronegativity
results[0] = electronegativity.calculateSigmaElectronegativity(container, atom);
// piElectronegativity
results[1] = pielectronegativity.calculatePiElectronegativity(container, atom);
// partialSigmaCharge
try {
peoe.assignGasteigerMarsiliSigmaPartialCharges(container, true);
} catch (Exception e) {
// ignored, underlying classes are logging this and this class
// is deprecated
}
results[2] = atom.getCharge();
// partialPiCharge
for (int i = 0; i < container.getAtomCount(); i++) container.getAtom(i).setCharge(0.0);
try {
pepe.assignGasteigerPiPartialCharges(container, true);
} catch (Exception e) {
LoggingToolFactory.createLoggingTool(IonizationPotentialTool.class).warn("Unexpected Error:", e);
}
results[3] = atom.getCharge();
// effectiveAtomicPolarizability
results[4] = pol.calculateGHEffectiveAtomPolarizability(container, atom, 100, true);
int position = container.indexOf(atom);
if (product != null)
results[5] = stabil.calculatePositive(product, product.getAtom(position));
else
results[5] = 0.0;
// numberResonance
IAtomContainer acR = gRI.getContainer(container, atom);
if (acR != null) {
results[6] = acR.getAtomCount();
// numberAromaticAtoms
// boolean isAromatic = Aromaticity.cdkLegacy().apply(container);
IRingSet ringSet = Cycles.sssr(container).toRingSet();
RingSetManipulator.markAromaticRings(ringSet);
int aromRingCount = 0;
for (IAtomContainer ring : ringSet.atomContainers()) {
if (ring.getFlag(CDKConstants.ISAROMATIC))
aromRingCount++;
}
results[7] = aromRingCount;
} else {
results[6] = 0;
results[7] = 0;
}
return results;
}
use of org.openscience.cdk.charges.Polarizability in project cdk by cdk.
the class IonizationPotentialTool method getQSARs.
/**
* Get the results of 7 qsar descriptors been applied. They are:
* Electronegativity,
* GasteigerMarsiliPartialCharges,
* GasteigerPEPEPartialCharges,
* Polarizability,
* StabilizationCharge,
* Number of Atom in resonance
* if the container in resonance is aromatic.
*
* @param container The IAtomContainer which contain the IAtom
* @param bond The IBond to calculate
* @return An Array containing the results
*/
public static double[] getQSARs(IAtomContainer container, IBond bond) throws CDKException {
Electronegativity electronegativity = new Electronegativity();
PiElectronegativity pielectronegativity = new PiElectronegativity();
GasteigerMarsiliPartialCharges peoe = new GasteigerMarsiliPartialCharges();
GasteigerPEPEPartialCharges pepe = new GasteigerPEPEPartialCharges();
Polarizability pol = new Polarizability();
StabilizationCharges stabil = new StabilizationCharges();
StructureResonanceGenerator gRI = new StructureResonanceGenerator();
double[] results = new double[7];
for (int ia = 0; ia < 2; ia++) {
IAtom atom = bond.getAtom(ia);
IAtomContainer product = initiateIonization(container, atom);
// sigmaElectronegativity
results[0] += electronegativity.calculateSigmaElectronegativity(container, atom);
// piElectronegativity
results[1] += pielectronegativity.calculatePiElectronegativity(container, atom);
// partialSigmaCharge
try {
peoe.assignGasteigerMarsiliSigmaPartialCharges(container, true);
} catch (Exception e) {
LoggingToolFactory.createLoggingTool(IonizationPotentialTool.class).warn("Unexpected Error:", e);
}
results[2] += atom.getCharge();
// partialPiCharge
for (int i = 0; i < container.getAtomCount(); i++) container.getAtom(i).setCharge(0.0);
try {
pepe.assignGasteigerPiPartialCharges(container, true);
} catch (Exception e) {
LoggingToolFactory.createLoggingTool(IonizationPotentialTool.class).warn("Unexpected Error:", e);
}
results[3] += atom.getCharge();
// effectiveAtomicPolarizability
results[4] += pol.calculateGHEffectiveAtomPolarizability(container, atom, 100, true);
int position = container.indexOf(atom);
if (product != null)
results[5] += stabil.calculatePositive(product, product.getAtom(position));
else
results[5] += 0.0;
// numberResonance
IAtomContainer acR = gRI.getContainer(container, atom);
if (acR != null) {
results[6] += acR.getAtomCount();
// numberAromaticAtoms
// boolean isAromatic = Aromaticity.cdkLegacy().apply(container);
// if(isAromatic)
// results[7] += 0.1;
} else {
results[6] += 0;
// results[7] += 0;
}
}
for (int i = 0; i < results.length; i++) results[i] = results[i] / 2;
return results;
}
use of org.openscience.cdk.charges.Polarizability in project cdk by cdk.
the class AutocorrelationDescriptorPolarizability method listpolarizability.
private static double[] listpolarizability(IAtomContainer container, int[][] dmat) throws CDKException {
int natom = container.getAtomCount();
double[] polars = new double[natom];
Polarizability polar = new Polarizability();
for (int i = 0; i < natom; i++) {
IAtom atom = container.getAtom(i);
try {
polars[i] = polar.calculateGHEffectiveAtomPolarizability(container, atom, false, dmat);
} catch (Exception ex1) {
throw new CDKException("Problems with assign Polarizability due to " + ex1, ex1);
}
}
return polars;
}
Aggregations