ipatomicdescriptor.java

化学图形处理软件

/* $Revision: 6228 $ $Author: egonw $ $Date: 2006-05-11 18:34:42 +0200 (Thu, 11 May 2006) $
 *
 * Copyright (C) 2006-2007  Miguel Rojas <[e-mail]>
 * 
 * Contact: [e-mail]
 *
 * This program is free software; you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation; either version 2 of the License, or
 * (at your option) any later version.
 *
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License along
 * with this program; if not, write to the Free Software Foundation, Inc.,
 * 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
 */
package org.openscience.cdk.qsar.descriptors.atomic;

import java.util.Iterator;
import java.util.List;

import org.openscience.cdk.AtomContainerSet;
import org.openscience.cdk.CDKConstants;
import org.openscience.cdk.aromaticity.HueckelAromaticityDetector;
import org.openscience.cdk.exception.CDKException;
import org.openscience.cdk.graph.invariant.ConjugatedPiSystemsDetector;
import org.openscience.cdk.interfaces.IAtom;
import org.openscience.cdk.interfaces.IAtomContainer;
import org.openscience.cdk.interfaces.IBond;
import org.openscience.cdk.interfaces.IMolecule;
import org.openscience.cdk.interfaces.IMoleculeSet;
import org.openscience.cdk.interfaces.IReaction;
import org.openscience.cdk.interfaces.IReactionSet;
import org.openscience.cdk.qsar.DescriptorSpecification;
import org.openscience.cdk.qsar.DescriptorValue;
import org.openscience.cdk.qsar.IAtomicDescriptor;
import org.openscience.cdk.qsar.descriptors.bond.BondPartialSigmaChargeDescriptor;
import org.openscience.cdk.qsar.descriptors.bond.ResonancePositiveChargeDescriptor;
import org.openscience.cdk.qsar.result.DoubleArrayResult;
import org.openscience.cdk.qsar.result.DoubleResult;
import org.openscience.cdk.reaction.IReactionProcess;
import org.openscience.cdk.reaction.type.ElectronImpactNBEReaction;

/**
 *  This class returns the ionization potential of an atom containg lone 
 *  pair electrons. It is
 *  based on a decision tree which is extracted from Weka(J48) from 
 *  experimental values. Up to now is only possible predict for 
 *  Cl,Br,I,N,P,O,S Atoms and they are not belong to
 *  conjugated system or not adjacent to an double bond.
 *
 * <p>This descriptor uses these parameters:
 * <table border="1">
 *   <tr>
 *     <td>Name</td>
 *     <td>Default</td>
 *     <td>Description</td>
 *   </tr>
 *   <tr>
 *     <td></td>
 *     <td></td>
 *     <td>no parameters</td>
 *   </tr>
 * </table>
 *
 * @author           Miguel Rojas
 * @cdk.created      2006-05-26
 * @cdk.module       qsar
 * @cdk.set          qsar-descriptors
 * @cdk.dictref      qsar-descriptors:ionizationPotential
 * 
 */
public class IPAtomicDescriptor implements IAtomicDescriptor {

	private IReactionSet reactionSet;
	
	/**
	 *  Constructor for the IPAtomicDescriptor object
	 */
	public IPAtomicDescriptor() {
	}
	/**
	 *  Gets the specification attribute of the IPAtomicDescriptor object
	 *
	 *@return    The specification value
	 */
	public DescriptorSpecification getSpecification() {
		return new DescriptorSpecification(
				"http://www.blueobelisk.org/ontologies/chemoinformatics-algorithms/#ionizationPotential",
				this.getClass().getName(),
				"$Id: IPAtomicDescriptor.java 6171 2006-5-22 19:29:58Z egonw $",
				"The Chemistry Development Kit");
	}
    /**
     * This descriptor does have any parameter.
     */
    public void setParameters(Object[] params) throws CDKException {
    }


    /**
     *  Gets the parameters attribute of the IPAtomicDescriptor object.
     *
     *@return    The parameters value
     * @see #setParameters
     */
    public Object[] getParameters() {
        return null;
    }
	/**
	 *  This method calculates the ionization potential of an atom.
	 *
	 *@param  chemObj           The IAtom to ionize.
	 *@param  container         Parameter is the IAtomContainer.
	 *@return                   The ionization potential. Not possible the ionization.
	 *@exception  CDKException  Description of the Exception
	 */
	public DescriptorValue calculate(IAtom atom, IAtomContainer atomContainer) throws CDKException{
        IAtomContainer localClone;
        try {
            localClone= (IAtomContainer) atomContainer.clone();
        } catch (CloneNotSupportedException e) {
            throw new CDKException("Error during clone");
        }
        IAtom clonedAtom = localClone.getAtom(atomContainer.getAtomNumber(atom));
        
        reactionSet = localClone.getBuilder().newReactionSet();
    	
		double resultD = -1.0;
		boolean isTarget = false;
		boolean isConjugated = false;
		double[] resultsH = null;

		/*control if it is into an aromatic or conjugated system*/
		HueckelAromaticityDetector.detectAromaticity(localClone,true);
		AtomContainerSet conjugatedPi = ConjugatedPiSystemsDetector.detect(localClone);
		Iterator acI = conjugatedPi.atomContainers();
 		while(acI.hasNext()){
			IAtomContainer ac = (IAtomContainer) acI.next();
			if(ac.contains(clonedAtom)){
				isConjugated = true;
				if(localClone.getMaximumBondOrder(clonedAtom) == 1 && localClone.getConnectedLonePairsCount(clonedAtom) > 0){
					resultsH = calculateHeteroAtomConjugatedDescriptor(clonedAtom, localClone,ac);
					resultD = getTreeHeteroConjAtom(resultsH);
					resultD += 0.05;
					isTarget = true;
				}
			}
 		}
 		
		if(atom.getFlag(CDKConstants.ISAROMATIC))
			return null;
		
		int count = localClone.getConnectedLonePairsCount(clonedAtom);
		
		if(localClone.getMaximumBondOrder(clonedAtom) > 1 && localClone.getConnectedLonePairsCount(clonedAtom) > 0){
			resultsH = calculateCarbonylDescriptor(clonedAtom, localClone);
			resultD = getTreeDoubleHetero(resultsH);
			resultD += 0.05;
			isTarget = true;
		}else if(localClone.getConnectedLonePairsCount(clonedAtom) > 0 && !isConjugated){
			resultsH = calculateHeteroAtomDescriptor(clonedAtom, localClone);
			resultD = getTreeHeteroAtom(resultsH);
			resultD += 0.05;
			isTarget = true;
		}
			
			
		if(isTarget){
			/* inizate reaction*/
			if(localClone.getConnectedLonePairsCount(clonedAtom) > 0){
				
				IMoleculeSet setOfReactants = localClone.getBuilder().newMoleculeSet();
				setOfReactants.addMolecule((IMolecule) localClone);
				IReactionProcess type  = new ElectronImpactNBEReaction();
				atom.setFlag(CDKConstants.REACTIVE_CENTER,true);
		        Object[] params = {Boolean.TRUE};
		        type.setParameters(params);
		        IReactionSet nbe = type.initiate(setOfReactants, null);
		        Iterator it = nbe.reactions();
		        while(it.hasNext()){
		        	IReaction reaction = (IReaction)it.next();
		        	reaction.setProperty("IonizationEnergy", new Double(resultD));
		        	reactionSet.addReaction(reaction);
		        }
			}
		}
		return new DescriptorValue(getSpecification(), getParameterNames(), getParameters(), new DoubleResult(resultD));
	}
	/**
	 * tree desicion for the carbonyl atoms
	 * 
	 * @param resultsH Array which contains the results of each descriptor
	 * @return the result
	 */
	private double getTreeHeteroConjAtom(double[] resultsH) {
		double result = 0.0;
		double SE_1 = resultsH[0];
		double SE_2 = resultsH[1];
		double EE_1  = resultsH[2];
		double RES_c2 = resultsH[3];
	
		if (SE_1 <= -0.069065)
		{
		  if (SE_1 <= -0.138994)
		  {
		    if (EE_1 <= 0.676304)
		    {
		      if (SE_1 <= -2.022267) { result = 07.7; /* 2.0 */}
		      else if (SE_1 > -2.022267) { result = 08.2; /* 4.0/2.0 */}
		    }
		    if (EE_1 > 0.676304)
		    {
		      if (SE_1 <= -2.514435) { result = 08.0; /* 3.0/1.0 */}
		      else if (SE_1 > -2.514435) { result = 09.5; /* 4.0/2.0 */}
		    }
		  }
		  if (SE_1 > -0.138994)
		  {
		    if (EE_1 <= 0.032157)
		    {
		      if (EE_1 <= 0.026958)
		      {
		        if (EE_1 <= 0.023022) { result = 08.1; /* 2.0/1.0 */}
		        else if (EE_1 > 0.023022)
		        {
		          if (SE_1 <= -0.074747) { result = 08.6; /* 3.0/1.0 */}
		          else if (SE_1 > -0.074747) { result = 08.9; /* 3.0 */}
		        }
		      }
		      if (EE_1 > 0.026958)
		      {
		        if (SE_2 <= 0.081676) { result = 08.3; /* 2.0/1.0 */}
		        else if (SE_2 > 0.081676) { result = 08.0; /* 6.0/1.0 */}
		      }
		    }
		    if (EE_1 > 0.032157)
		    {
		      if (EE_1 <= 0.045671) { result = 08.6; /* 5.0/2.0 */}
		      else if (EE_1 > 0.045671)
		      {
		        if (SE_1 <= -0.137681) { result = 08.1; /* 2.0/1.0 */}
		        else if (SE_1 > -0.137681) { result = 08.4; /* 5.0/2.0 */}
		      }
		    }
		  }
		}
		if (SE_1 > -0.069065)
		{
		  if (SE_1 <= -0.061906)
		  {
		    if (RES_c2 <= 0.006978) { result = 07.4; /* 15.0/10.0 */}
		    else if (RES_c2 > 0.006978)
		    {
		      if (SE_1 <= -0.063064) { result = 07.5; /* 8.0/5.0 */}
		      else if (SE_1 > -0.063064) { result = 07.3; /* 2.0/1.0 */}
		    }
		  }
		  if (SE_1 > -0.061906)
		  {
		    if (RES_c2 <= 0.003398)
		    {
		      if (SE_1 <= 0.051565)
		      {
		        if (EE_1 <= 0.007932)
		        {
		          if (SE_1 <= -0.023796) { result = 08.7; /* 4.0/2.0 */}
		          else if (SE_1 > -0.023796)
		          {
		            if (EE_1 <= -0.003129) { result = 07.3; /* 2.0/1.0 */}
		            else if (EE_1 > -0.003129) { result = 07.4; /* 4.0/2.0 */}
		          }
		        }
		        if (EE_1 > 0.007932)
		        {
		          if (SE_1 <= -0.0252)
		          {
		            if (SE_1 <= -0.02952) { result = 08.0; /* 2.0/1.0 */}
		            else if (SE_1 > -0.02952) { result = 07.8; /* 2.0/1.0 */}
		          }
		          if (SE_1 > -0.0252) { result = 07.7; /* 4.0/1.0 */}
		        }
		      }
		      if (SE_1 > 0.051565)
		      {
		        if (RES_c2 <= -0.003572) { result = 08.5; /* 4.0/2.0 */}
		        else if (RES_c2 > -0.003572)
		        {
		          if (SE_1 <= 0.126992) { result = 07.6; /* 2.0/1.0 */}
		          else if (SE_1 > 0.126992) { result = 08.2; /* 2.0/1.0 */}
		        }
		      }
		    }
		    if (RES_c2 > 0.003398)
		    {
		      if (SE_1 <= -0.031039)
		      {
		        if (RES_c2 <= 0.005076) { result = 08.1; /* 9.0/5.0 */}
		        else if (RES_c2 > 0.005076)
		        {
		          if (SE_1 <= -0.061705) { result = 08.2; /* 2.0 */}
		          else if (SE_1 > -0.061705)
		          {
		            if (SE_1 <= -0.060998) { result = 08.3; /* 2.0 */}
		            else if (SE_1 > -0.060998) { result = 08.0; /* 2.0/1.0 */}
		          }
		        }
		      }
		      if (SE_1 > -0.031039)
		      {
		        if (SE_1 <= -0.023566)
		        {
		          if (SE_1 <= -0.028788) { result = 07.9; /* 2.0 */}
		          else if (SE_1 > -0.028788) { result = 09.3; /* 2.0/1.0 */}
		        }
		        if (SE_1 > -0.023566)
		        {
		          if (RES_c2 <= 0.004518) { result = 08.2; /* 2.0/1.0 */}
		          else if (RES_c2 > 0.004518) { result = 08.7; /* 4.0/2.0 */}
		        }
		      }
		    }
		  }
		}
		
		return result;
	}
	/**
	 * tree desicion for the carbonyl atoms
	 * 
	 * @param resultsH Array which contains the results of each descriptor
	 * @return the result
	 */
	private double getTreeDoubleHetero(double[] resultsH) {
		double result = 0.0;
		double SE_c = resultsH[0];
		double PCH_c = resultsH[1];
		double SB  = resultsH[2];
		double SE_x = resultsH[3];
		double PCH_x = resultsH[4];
		double RES_c = resultsH[5];
		
		if (PCH_c <= 0.045111)
		{
		  if (PCH_x <= -0.041368)
		  {
		    if (SE_c <= 7.471265)
		    {
		      if (RES_c <= 0.182146)
		      {
		        if (SB <= 0.448987)
		        {
		          if (RES_c <= 0.181878) { result = 08.4; /* 3.0/1.0 */}
		          else if (RES_c > 0.181878) { result = 09.2; /* 4.0/1.0 */}
		        }
		        if (SB > 0.448987)
		        {
		          if (SE_c <= 6.758953) { result = 09.0; /* 2.0 */}
		          else if (SE_c > 6.758953) { result = 09.6; /* 2.0/1.0 */}
		        }
		      }
		      if (RES_c > 0.182146)
		      {
		        if (RES_c <= 0.371342)
		        {
		          if (PCH_c <= 0.019491) { result = 09.3; /* 9.0/3.0 */}
		          else if (PCH_c > 0.019491)
		          {
		            if (SE_c <= 6.683097) { result = 09.7; /* 3.0 */}
		            else if (SE_c > 6.683097) { result = 09.6; /* 4.0/2.0 */}
		          }
		        }
		        if (RES_c > 0.371342)
		        {
		          if (PCH_c <= 0.019061) { result = 08.7; /* 3.0/1.0 */}
		          else if (PCH_c > 0.019061)
		          {
		            if (SE_c <= 6.689165) { result = 09.9; /* 2.0/1.0 */}
		            else if (SE_c > 6.689165) { result = 09.5; /* 4.0/1.0 */}
		          }
		        }
		      }
		    }
		    if (SE_c > 7.471265)
		    {
		      if (SE_x <= 13.191725)
		      {
		        if (SB <= 0.479515)
		        {
		          if (PCH_x <= -0.045111) { result = 10.2; /* 4.0/2.0 */}
		          else if (PCH_x > -0.045111) { result = 09.8; /* 2.0/1.0 */}
		        }
		        if (SB > 0.479515)
		        {
		          if (SE_c <= 7.571876) { result = 10.0; /* 4.0/2.0 */}