gasteigerpepepartialcharges.java

化学图形处理软件

    			mol.getBond(k).setFlag(CDKConstants.REACTIVE_CENTER,false);
    			mol.getBond(k).getAtom(0).setFlag(CDKConstants.REACTIVE_CENTER,false);
    			mol.getBond(k).getAtom(1).setFlag(CDKConstants.REACTIVE_CENTER,false);
    		}
    		setOfM2.addMolecule((IMolecule) mol);
    		Object[] params2 = {Boolean.FALSE};
			type.setParameters(params2);
			IReactionSet setOfReactions2 = type.initiate(setOfM2, null);
			if(setOfReactions2.getReactionCount() > 0){
				
			IMolecule react = setOfReactions2.getReaction(0).getReactants().getMolecule(0);

			set.addAtomContainer(react);
			}
        }

		return set;
	}
	/**
	 * get the electrostatic potential of the neighbours of a atom.
	 *  
	 * @param ac   The IAtomContainer to study
	 * @param ds 
	 * @param atom1 The position of the IAtom to study
	 * @return     The sum of electrostatic potential of the neighbours
	 */
	private double getElectrostaticPotentialN(IAtomContainer ac, int atom1, double[] ds) {
		
//		double CoulombForceConstant = 1/(4*Math.PI*8.81/*Math.pow(10, -12)*/);
		double CoulombForceConstant = 0.048;
		double sum = 0.0;
		try {
			if (factory == null) 
                factory = AtomTypeFactory.getInstance(
                    "org/openscience/cdk/config/data/jmol_atomtypes.txt", 
                    ac.getBuilder()
                );

			java.util.List atoms = ac.getConnectedAtomsList(ac.getAtom(atom1));
			for(int i = 0 ; i < atoms.size() ; i++){
				double covalentradius = 0;
	            String symbol = ((IAtom)atoms.get(i)).getSymbol();
	            IAtomType type = factory.getAtomType(symbol);
	            covalentradius = type.getCovalentRadius();

	            double charge = ds[STEP_SIZE * atom1 + atom1 + 5];
				double sumI = CoulombForceConstant*charge/(covalentradius*covalentradius);
//				logger.debug("sum_("+sumI+") = CFC("+CoulombForceConstant+")*charge("+charge+"/ret("+covalentradius);
				sum += sumI;
			}
		} catch (CDKException e) {
            logger.debug(e);
        }
		
		return sum;
	}


	/**
	 * get the topological weight factor for each atomContainer
	 * 
	 * @param atomContainer  The IAtomContainer to study.
	 * @param ac             The IAtomContainer to study.
	 * @return The value
	 */
	private double getTopologicalFactors(IAtomContainer atomContainer,IAtomContainer ac) {
		/*factor for separation of charge*/
		int totalNCharge1 = AtomContainerManipulator.getTotalNegativeFormalCharge(atomContainer);
		int totalPCharge1 = AtomContainerManipulator.getTotalPositiveFormalCharge(atomContainer);
		
		double fQ = 1.0;
		if(totalNCharge1 != 0.0){
			fQ = 0.5;
			for(int i = 0; i < atomContainer.getBondCount(); i++){
				IBond bond = atomContainer.getBond(i);
				if(bond.getAtom(0).getFormalCharge() != 0.0 && bond.getAtom(1).getFormalCharge() != 0.0){
					fQ = 0.25;
					break;
				}
			}
		}
		/*factor, if the number of covalents bonds is decreased*/
		double fB = 1.0;
		
		int numBond1 = 0;
		int numBond2 = 0;
        for (int i = 0; i < atomContainer.getBondCount(); i++) {
            if (atomContainer.getBond(i).getOrder() == 2.0)
                numBond1 += 1;
            if (ac.getBond(i).getOrder() == 2.0)
                numBond2 += 1;
        }
        
        if(numBond1 </*>*/ numBond2)
			fB = 0.8;
        
		double fPlus = 1.0;
        if(totalNCharge1 == 0.0 && totalPCharge1 == 0.0 )
        	fPlus = 0.1;
        
        
        /*aromatic*/
        double fA = 1.0;
        try {
			if(HueckelAromaticityDetector.detectAromaticity(ac))
				if(!HueckelAromaticityDetector.detectAromaticity(atomContainer))
						fA = 0.3;
		} catch (CDKException e) {
			e.printStackTrace();
		}
//		logger.debug("return "+fQ*fB*fPlus*fA+"= sp:"+fQ+", dc:"+fB+", fPlus:"+fPlus+", fA:"+fA);
		
		return fQ*fB*fPlus*fA;
	}


	/**
	 *  Get the StepSize attribute of the GasteigerMarsiliPartialCharges
	 *  object
	 *
	 *@return STEP_SIZE
	 */
	public int getStepSize(){
		return STEP_SIZE;
	}
	

	/**
	 * Method which stores and assigns the factors a,b,c and CHI+
	 *
	 * @return     Array of doubles [a1,b1,c1,denom1,chi1,q1...an,bn,cn...] 1:Atom 1-n in AtomContainer
	 */
	private double[][] assignPiFactors(IAtomContainerSet setAc) {
		//a,b,c,denom,chi,q
		double[][] gasteigerFactors = new double[setAc.getAtomContainerCount()][(setAc.getAtomContainer(0).getAtomCount() * (STEP_SIZE+1))];
		String AtomSymbol = "";
		double[] factors = new double[]{0.0, 0.0, 0.0};
		for( int k = 1 ; k < setAc.getAtomContainerCount(); k ++){
			IAtomContainer ac = setAc.getAtomContainer(k);
			for (int i = 0; i < ac.getAtomCount(); i++) {
				factors[0] = 0.0;
				factors[1] = 0.0;
				factors[2] = 0.0;
				AtomSymbol = ac.getAtom(i).getSymbol();
				if (AtomSymbol.equals("H")) {
					factors[0] = 0.0;
					factors[1] = 0.0;
					factors[2] = 0.0;
				} else if (AtomSymbol.equals("C")) {/*
					if(ac.getAtom(i).getFlag(ISCHANGEDFC))*/{
						factors[0] = 5.60;
						factors[1] = 8.93;
						factors[2] = 2.94;
					}
				} else if (AtomSymbol.equals("O")) {
						if(ac.getMaximumBondOrder(ac.getAtom(i)) == 1){
							factors[0] = 10.0;
							factors[1] = 13.86;
							factors[2] = 9.68;
						}else {
							factors[0] = 7.91;
							factors[1] = 14.76;
							factors[2] = 6.85;
						}
				} else if (AtomSymbol.equals("N")) {
					if(ac.getMaximumBondOrder(ac.getAtom(i)) != 1){
						factors[0] = 7.95;/*7.95*/
						factors[1] = 9.73;/*9.73*/
						factors[2] = 2.67;/*2.67*/
					}else {
						factors[0] = 4.54;/*4.54*//*5.5*/
						factors[1] = 11.86;/*11.86*//*10.86*/
						factors[2] = 7.32;/*7.32*//*7.99*/
					}
				} else if (AtomSymbol.equals("S")) {
					if(ac.getMaximumBondOrder(ac.getAtom(i)) == 1){
						factors[0] = 7.73;
						factors[1] = 8.16;
						factors[2] = 1.81;
					}else {
						factors[0] = 6.60;
						factors[1] = 10.32;
						factors[2] = 3.72;
					}
				} else if (AtomSymbol.equals("F")) {
					factors[0] = 7.34;
					factors[1] = 13.86;
					factors[2] = 9.68;
				} else if (AtomSymbol.equals("Cl")) {
					factors[0] = 6.50;
					factors[1] = 11.02;
					factors[2] = 4.52;
				} else if (AtomSymbol.equals("Br")) {
					factors[0] = 5.20;
					factors[1] = 9.68;
					factors[2] = 4.48;
				} else if (AtomSymbol.equals("I")) {
					factors[0] = 4.95;
					factors[1] = 8.81;
					factors[2] = 3.86;
				}
			
				gasteigerFactors[k][STEP_SIZE * i + i] = factors[0];
				gasteigerFactors[k][STEP_SIZE * i + i + 1] = factors[1];
				gasteigerFactors[k][STEP_SIZE * i + i + 2] = factors[2];
				gasteigerFactors[k][STEP_SIZE * i + i + 5] = ac.getAtom(i).getCharge();
	
				if (factors[0] == 0 && factors[1] == 0 && factors[2] == 0) {
					gasteigerFactors[k][STEP_SIZE * i + i + 3] = 1;
				} else {
					gasteigerFactors[k][STEP_SIZE * i + i + 3] = factors[0] + factors[1] + factors[2];
				}
			}
		}
		

		return gasteigerFactors;
	}
	/**
	 *  Method which stores and assigns the factors a,b,c and CHI+
	 *
	 *@return     Array of doubles [a1,b1,c1,denom1,chi1,q1...an,bn,cn...] 1:Atom 1-n in AtomContainer
	 */
	public double[][] assignrPiMarsilliFactors(IAtomContainerSet setAc) {
		//a,b,c,denom,chi,q
		double[][] gasteigerFactors = new double[setAc.getAtomContainerCount()][(setAc.getAtomContainer(0).getAtomCount() * (STEP_SIZE+1))];
		String AtomSymbol = "";
		double[] factors = new double[]{0.0, 0.0, 0.0};
		for( int k = 1 ; k < setAc.getAtomContainerCount(); k ++){
			IAtomContainer ac = setAc.getAtomContainer(k);
			
			for (int i = 0; i < ac.getAtomCount(); i++) {
				factors[0] = 0.0;
				factors[1] = 0.0;
				factors[2] = 0.0;
				AtomSymbol = ac.getAtom(i).getSymbol();
				if (AtomSymbol.equals("H")) {
					factors[0] = 0.0;
					factors[1] = 0.0;
					factors[2] = 0.0;
				} else if (AtomSymbol.equals("C")) {
						factors[0] = 5.98;/*5.98-5.60*/
						factors[1] = 7.93;/*7.93-8.93*/
						factors[2] = 1.94;
				} else if (AtomSymbol.equals("O")) {
						if(ac.getMaximumBondOrder(ac.getAtom(i)) > 1){
							factors[0] = 11.2;/*11.2-10.0*/
							factors[1] = 13.24;/*13.24-13.86*/
							factors[2] = 9.68;
						}else {
							factors[0] = 7.91;
							factors[1] = 14.76;
							factors[2] = 6.85;
						}
				} else if (AtomSymbol.equals("N")) {
					if(ac.getMaximumBondOrder(ac.getAtom(i))  > 1){
						
						factors[0] = 8.95;/*7.95*/
						factors[1] = 9.73;/*9.73*/
						factors[2] = 2.67;/*2.67*/
					}else {
						factors[0] = 4.54;
						factors[1] = 11.86;
						factors[2] = 7.32;
					}
				} else if (AtomSymbol.equals("P")) {// <--No correct
					if(ac.getMaximumBondOrder(ac.getAtom(i))  > 1){
						factors[0] = 10.73;// <--No correct
						factors[1] = 11.16;// <--No correct
						factors[2] = 6.81;// <--No correct
					}else {
						factors[0] = 9.60;// <--No correct
						factors[1] = 13.32;// <--No correct
						factors[2] = 2.72;// <--No correct
					}
				} else if (AtomSymbol.equals("S")) {
					if(ac.getMaximumBondOrder(ac.getAtom(i))  > 1){
						
						factors[0] = 7.73;
						factors[1] = 8.16;
						factors[2] = 1.81;
					}else {
						factors[0] = 6.60;
						factors[1] = 10.32;
						factors[2] = 3.72;
					}
				} else if (AtomSymbol.equals("F")) {
					factors[0] = 7.14/*7.34*/;
					factors[1] = 13.86;
					factors[2] = 5.68;
				} else if (AtomSymbol.equals("Cl")) {
					factors[0] = 6.51;/*6.50*/
					factors[1] = 11.02;
					factors[2] = 4.52;
				} else if (AtomSymbol.equals("Br")) {
					factors[0] = 5.20;
					factors[1] = 9.68;
					factors[2] = 4.48;
				} else if (AtomSymbol.equals("I")) {
					factors[0] = 4.95;
					factors[1] = 8.81;
					factors[2] = 3.86;
				}
			
				gasteigerFactors[k][STEP_SIZE * i + i] = factors[0];
				gasteigerFactors[k][STEP_SIZE * i + i + 1] = factors[1];
				gasteigerFactors[k][STEP_SIZE * i + i + 2] = factors[2];
				gasteigerFactors[k][STEP_SIZE * i + i + 5] = ac.getAtom(i).getCharge();
				
				if (factors[0] == 0 && factors[1] == 0 && factors[2] == 0) {
					gasteigerFactors[k][STEP_SIZE * i + i + 3] = 1;
				} else {
					gasteigerFactors[k][STEP_SIZE * i + i + 3] = factors[0] + factors[1] + factors[2];
				}
			}
		}
		

		return gasteigerFactors;
	}
	/**
     * clean the flags CDKConstants.REACTIVE_CENTER from the molecule
     * 
	 * @param ac
	 */
	public void cleanFlagReactiveCenter(IAtomContainer ac){
		for(int j = 0 ; j < ac.getAtomCount(); j++)
			ac.getAtom(j).setFlag(CDKConstants.REACTIVE_CENTER, false);
		for(int j = 0 ; j < ac.getBondCount(); j++)
			ac.getBond(j).setFlag(CDKConstants.REACTIVE_CENTER, false);
	}
}