PDBlib Class Reference

PDBlib Class. More...

#include <PDBlib.hh>

Collaboration diagram for PDBlib:

Public Member Functions
	PDBlib ()
	First constructor. More...
	~PDBlib ()
	Destructor. More...
Molecule *	Load (const string &filename, unsigned short int &isDNA, unsigned short int verbose)
	Load PDB file into memory. More...
Barycenter *	ComputeNucleotideBarycenters (Molecule *moleculeListTemp)
	Compute nucleotide barycenter from memory. More...
void	ComputeBoundingVolumeParams (Molecule *moleculeListTemp, double &dX, double &dY, double &dZ, double &tX, double &tY, double &tZ)
	Compute the corresponding bounding volume parameters. More...
void	ComputeNbNucleotidsPerStrand (Molecule *moleculeListTemp)
	Compute number of nucleotide per strand. More...
unsigned short int	ComputeMatchEdepDNA (Barycenter *, Molecule *, double x, double y, double z, int &numStrand, int &numNucleotid, int &codeResidue)
	Compute if energy is deposited in per atom. More...

Private Member Functions
double	DistanceTwo3Dpoints (double xA, double xB, double yA, double yB, double zA, double zB)
	return distance between two 3D points More...

Private Attributes
int	fNbNucleotidsPerStrand
	Number of nucleotid per strand. More...

Detailed Description

PDBlib Class.

This Class define Molecule model ...

Definition at line 58 of file PDBlib.hh.

Constructor & Destructor Documentation

PDBlib()

PDBlib::PDBlib

(

)

First constructor.

Definition at line 58 of file PDBlib.cc.

              :fNbNucleotidsPerStrand(0)
{
}

~PDBlib()

PDBlib::~PDBlib ( )

inline

Destructor.

Definition at line 64 of file PDBlib.hh.

{};

Here is the call graph for this function:

Member Function Documentation

ComputeBoundingVolumeParams()

void PDBlib::ComputeBoundingVolumeParams	(	Molecule *	moleculeListTemp,
		double &	dX,
		double &	dY,
		double &	dZ,
		double &	tX,
		double &	tY,
		double &	tZ
	)

Compute the corresponding bounding volume parameters.

the corresponding bounding volume parameters

the corresponding bounding volume parameters to build a box from atoms coordinates

Definition at line 667 of file PDBlib.cc.

{
  double minminX,minminY,minminZ;    //minimum minimorum
  double maxmaxX,maxmaxY,maxmaxZ;    //maximum maximorum

#ifdef GEANT4
  minminX=2000000000;
  minminY=2000000000;
  minminZ=2000000000;
  maxmaxX=-2000000000;
  maxmaxY=-2000000000;
  maxmaxZ=-2000000000;
#else
  minminX=numeric_limits<double>::max();
  minminY=numeric_limits<double>::max();
  minminZ=numeric_limits<double>::max();
  maxmaxX=numeric_limits<double>::min();
  maxmaxY=numeric_limits<double>::min();
  maxmaxZ=numeric_limits<double>::min();
#endif

  while (moleculeListTemp)
  {
    if (minminX > moleculeListTemp->fMaxGlobX)
    {
      minminX = moleculeListTemp->fMaxGlobX;
    }
    if (minminY > moleculeListTemp->fMaxGlobY)
    {
      minminY = moleculeListTemp->fMaxGlobY;
    }
    if (minminZ > moleculeListTemp->fMaxGlobZ)
    {
      minminZ = moleculeListTemp->fMaxGlobZ;
    }
    if (maxmaxX < moleculeListTemp->fMinGlobX)
    {
      maxmaxX = moleculeListTemp->fMinGlobX;
    }
    if (maxmaxY < moleculeListTemp->fMinGlobY)
    {
      maxmaxY = moleculeListTemp->fMinGlobY;
    }
    if (maxmaxZ < moleculeListTemp->fMinGlobZ)
    {
      maxmaxZ = moleculeListTemp->fMinGlobZ;
    }

    moleculeListTemp=moleculeListTemp->GetNext();
  }//end of while sur moleculeListTemp

  dX=(maxmaxX-minminX)/2.+1.8;  //1.8 => size of biggest radius for atoms
  dY=(maxmaxY-minminY)/2.+1.8;
  dZ=(maxmaxZ-minminZ)/2.+1.8;

  tX=minminX+(maxmaxX-minminX)/2.;
  tY=minminY+(maxmaxY-minminY)/2.;
  tZ=minminZ+(maxmaxZ-minminZ)/2.;
}

Here is the call graph for this function:

Here is the caller graph for this function:

ComputeMatchEdepDNA()

unsigned short int PDBlib::ComputeMatchEdepDNA	(	Barycenter *	BarycenterList,
		Molecule *	moleculeListTemp,
		double	x,
		double	y,
		double	z,
		int &	numStrand,
		int &	numNucleotid,
		int &	codeResidue
	)

Compute if energy is deposited in per atom.

Compute barycenters.

Compute barycenters and its coordinate for nucleotides

Parameters

Molecule

* moleculeList

Returns

Barycenter *

Definition at line 772 of file PDBlib.cc.

{
  unsigned short int matchFound=0;

  Molecule *mLTsavedPointer = moleculeListTemp;
  Barycenter *BLsavedPointer = BarycenterList;

  short int strandNum=0;//Strand number
  int residueNum=1;//Residue (nucleotide) number
  string baseName;//Base name [A,C,T,G]
  unsigned short int BSP=2;//Base (default value), Sugar, Phosphat

  double smallestDist;//smallest dist Atom <-> edep coordinates
  double distEdepDNA;
  double distEdepAtom;

  //Residue (Base, Phosphate,suggar) list
  Residue *residueListTemp;
  //Atom list inside a residue
  Atom *AtomTemp;

  int k = 0;//Molecule number
  moleculeListTemp = mLTsavedPointer;
  BarycenterList = BLsavedPointer;

  smallestDist=33.0;//Sufficiently large value
  while (moleculeListTemp)
  {
    k++;
    residueListTemp = moleculeListTemp->GetFirst();

    int j = 0;//Residue number

#ifdef GEANT4
    int j_save=2000000000;//Saved res. number if match found
#else
    int j_save = numeric_limits<int>::max();//Saved res. number if match found
#endif

    while (residueListTemp)
    {
      j++;

      if (j - j_save > 2 ) break;

      distEdepDNA=DistanceTwo3Dpoints(x,BarycenterList->fCenterX,
                                      y,BarycenterList->fCenterY,
                                      z,BarycenterList->fCenterZ);
      if (distEdepDNA < BarycenterList->GetRadius())
      {
        //Find the closest atom
        //Compute distance between energy deposited and atoms for a residue
        //if distance <1.8 then match OK but search inside 2 next residues
        AtomTemp=residueListTemp->GetFirst();
        for (int iii=0 ; iii < residueListTemp->fNbAtom ; iii++)
        {

          distEdepAtom=DistanceTwo3Dpoints(x,AtomTemp->GetX(),
                                           y,AtomTemp->GetY(),
                                           z,AtomTemp->GetZ());

          if ((distEdepAtom < AtomTemp->GetVanDerWaalsRadius())
              && (smallestDist > distEdepAtom))
          {
            strandNum=k;

            if (k==2){
              residueNum = fNbNucleotidsPerStrand + 1 -
                  residueListTemp->fResSeq;
            }
            else
            {
              residueNum = residueListTemp->fResSeq;
            }

            baseName = (residueListTemp->fResName)[2];
            if (residueListTemp->fResSeq == 1)
            {
              if (iii == 0) BSP = 0;//"Phosphate"
              else if (iii < 8) BSP = 1;//"Sugar"
              else BSP = 2;//"Base"
            }
            else
            {
              if (iii < 4) BSP = 0;//"Phosphate"
              else if (iii < 11) BSP = 1;//"Sugar"
              else BSP = 2;//"Base"
            }

            smallestDist=distEdepAtom;

#ifdef GEANT4
            int j_max_value=2000000000;
#else
            int j_max_value = numeric_limits<int>::max();
#endif
            if (j_save == j_max_value) j_save = j;
            matchFound = 1;
          }
          AtomTemp=AtomTemp->GetNext();
        }//end of for (  iii=0 ; iii < residueListTemp->nbAtom ; iii++)
      }//end for if (distEdepDNA < BarycenterList->GetRadius())
      BarycenterList=BarycenterList->GetNext();
      residueListTemp=residueListTemp->GetNext();
    }//end of while sur residueListTemp
    moleculeListTemp=moleculeListTemp->GetNext();
  }//end of while sur moleculeListTemp

  numStrand = strandNum;
  numNucleotid = residueNum;
  codeResidue = BSP;

  return matchFound;
}

Here is the call graph for this function:

Here is the caller graph for this function:

ComputeNbNucleotidsPerStrand()

void PDBlib::ComputeNbNucleotidsPerStrand

(

Molecule *

moleculeListTemp

)

Compute number of nucleotide per strand.

Compute number of nucleotide per strand for DNA

Definition at line 736 of file PDBlib.cc.

{
  Residue *residueListTemp;

  int k=0;
  int j_old=0;

  while (moleculeListTemp)
  {
    residueListTemp = moleculeListTemp->GetFirst();

    k++;
    int j=0;

    while (residueListTemp)
    {
      j++;
      residueListTemp=residueListTemp->GetNext();
    }//end of while sur residueListTemp

    j_old+=j;
    moleculeListTemp=moleculeListTemp->GetNext();
  }//end of while sur moleculeListTemp

  fNbNucleotidsPerStrand=j_old/2;
}

Here is the call graph for this function:

Here is the caller graph for this function:

ComputeNucleotideBarycenters()

Barycenter * PDBlib::ComputeNucleotideBarycenters

(

Molecule *

moleculeListTemp

)

Compute nucleotide barycenter from memory.

Compute barycenters.

Compute barycenters and its coordinate for nucleotides

Parameters

Molecule

* moleculeList

Returns

Barycenter *

molecs->push_back(*moleculeListTemp);

Definition at line 480 of file PDBlib.cc.

{
  //Placement and physical volume construction from memory
  Barycenter * BarycenterFirst = NULL;
  Barycenter * BarycenterOld = NULL;
  Barycenter * BarycenterNext = NULL;

  //Residue (Base, Phosphate,sugar) list
  Residue *residueListTemp;
  //Atom list inside a residu
  Atom *AtomTemp;

  int k=0;
  int j_old=0;

  while (moleculeListTemp)
  {
    residueListTemp = moleculeListTemp->GetFirst();

    k++;
    int j=0;

    //Check numerotation style (1->n per strand or 1->2n for two strand)
    int correctNumerotationNumber=0;
    if (k==2 && residueListTemp->fResSeq > 1)
    {
      correctNumerotationNumber=residueListTemp->fResSeq;
    }

    while (residueListTemp)
    {
      AtomTemp=residueListTemp->GetFirst();
      j++;

      //Correction consequently to numerotation check
      if (correctNumerotationNumber!=0)
      {
        residueListTemp->fResSeq=residueListTemp->fResSeq-
            correctNumerotationNumber+1;
      }

      //Barycenter computation
      double baryX=0.,baryY=0.,baryZ=0.;
      double baryBaseX=0.,baryBaseY=0.,baryBaseZ=0.;
      double barySugX=0.,barySugY=0.,barySugZ=0.;
      double baryPhosX=0.,baryPhosY=0.,baryPhosZ=0.;
      unsigned short int nbAtomInBase=0;

      for (int i=0 ; i < residueListTemp->fNbAtom ; i++)
      {
        //Compute barycenter of the nucletotide
        baryX+=AtomTemp->fX;
        baryY+=AtomTemp->fY;
        baryZ+=AtomTemp->fZ;
        //Compute barycenters for Base Sugar Phosphat
        if (residueListTemp->fResSeq == 1)
        {
          if (i==0)
          {
            baryPhosX+=AtomTemp->fX;
            baryPhosY+=AtomTemp->fY;
            baryPhosZ+=AtomTemp->fZ;
          }
          else if (i<8)
          {
            barySugX+=AtomTemp->fX;
            barySugY+=AtomTemp->fY;
            barySugZ+=AtomTemp->fZ;
          }
          else
          {
            //hydrogen are placed at the end of the residue in a PDB file
            //We don't want them for this calculation
            if (AtomTemp->fElement!="H"){
              baryBaseX+=AtomTemp->fX;
              baryBaseY+=AtomTemp->fY;
              baryBaseZ+=AtomTemp->fZ;
              nbAtomInBase++;}
          }
        }
        else
        {
          if (i<4)
          {
            baryPhosX+=AtomTemp->fX;
            baryPhosY+=AtomTemp->fY;
            baryPhosZ+=AtomTemp->fZ;
          }
          else if (i<11)
          {
            barySugX+=AtomTemp->fX;
            barySugY+=AtomTemp->fY;
            barySugZ+=AtomTemp->fZ;
          }
          else
          {
            //hydrogen are placed at the end of the residue in a PDB file
            //We don't want them for this calculation
            if (AtomTemp->fElement!="H"){// break;
              baryBaseX+=AtomTemp->fX;
              baryBaseY+=AtomTemp->fY;
              baryBaseZ+=AtomTemp->fZ;
              nbAtomInBase++;}
          }
        }
        AtomTemp=AtomTemp->GetNext();
      }//end of for (  i=0 ; i < residueListTemp->nbAtom ; i++)

      baryX = baryX / (double)residueListTemp->fNbAtom;
      baryY = baryY / (double)residueListTemp->fNbAtom;
      baryZ = baryZ / (double)residueListTemp->fNbAtom;

      if (residueListTemp->fResSeq != 1) //Special case first Phosphate
      {
        baryPhosX = baryPhosX / 4.;
        baryPhosY = baryPhosY / 4.;
        baryPhosZ = baryPhosZ / 4.;
      }
      barySugX = barySugX / 7.;
      barySugY = barySugY / 7.;
      barySugZ = barySugZ / 7.;
      baryBaseX = baryBaseX / (double)nbAtomInBase;
      baryBaseY = baryBaseY / (double)nbAtomInBase;
      baryBaseZ = baryBaseZ / (double)nbAtomInBase;

      //Barycenter creation:
      if (BarycenterOld == NULL)
      {
        BarycenterOld =new Barycenter(j+j_old,baryX,baryY,baryZ, //j [1..n]
                                      baryBaseX,baryBaseY,baryBaseZ,
                                      barySugX,barySugY,barySugZ,
                                      baryPhosX,baryPhosY,baryPhosZ);
        BarycenterFirst = BarycenterOld;
      }
      else
      {
        BarycenterNext = new Barycenter(j+j_old,baryX,baryY,baryZ,
                                        baryBaseX,baryBaseY,baryBaseZ,
                                        barySugX,barySugY,barySugZ,
                                        baryPhosX,baryPhosY,baryPhosZ);
        BarycenterOld -> SetNext(BarycenterNext);
        BarycenterOld = BarycenterNext;
      }

      //distance computation between all atoms inside
      //a residue and the barycenter
      AtomTemp=residueListTemp->GetFirst();
      double dT3Dp;
      double max=0.;
      for (int ii=0 ; ii < residueListTemp->fNbAtom ; ii++)
      {
        dT3Dp = DistanceTwo3Dpoints(AtomTemp->fX,BarycenterOld->fCenterX,
                                    AtomTemp->fY,BarycenterOld->fCenterY,
                                    AtomTemp->fZ,BarycenterOld->fCenterZ);
        BarycenterOld->SetDistance(ii,dT3Dp);
        if (dT3Dp>max) max=dT3Dp;
        AtomTemp=AtomTemp->GetNext();
      }//end of for (  i=0 ; i < residueListTemp->nbAtom ; i++)

      BarycenterOld->SetRadius(max+1.8);
      residueListTemp=residueListTemp->GetNext();

    }//end of while sur residueListTemp

    j_old+=j;

    moleculeListTemp=moleculeListTemp->GetNext();
  }//end of while sur moleculeListTemp

  if(BarycenterNext!=NULL)
  {
    BarycenterNext -> SetNext(NULL);
  }

  return BarycenterFirst;
}

Here is the call graph for this function:

Here is the caller graph for this function:

DistanceTwo3Dpoints()

double PDBlib::DistanceTwo3Dpoints ( double  xA, double  xB, double  yA, double  yB, double  zA, double  zB  )

private

return distance between two 3D points

Definition at line 892 of file PDBlib.cc.

{
  return std::sqrt ( (xA-xB)*(xA-xB) + (yA-yB)*(yA-yB) + (zA-zB)*(zA-zB) );
}

Here is the caller graph for this function:

Load()

Molecule * PDBlib::Load	(	const string &	filename,
		unsigned short int &	isDNA,
		unsigned short int	verbose = 0
	)

Load PDB file into memory.

Load a PDB file into memory.

molecule (polymer,?), read line, key words

Parameters

filename	string for filename
verbosity

Returns

List of molecules

Definition at line 73 of file PDBlib.cc.

{
  string sLine = "";
  ifstream infile;
  infile.open(filename.c_str());
  if (!infile)
  {
#ifdef GEANT4
    G4cout<<"PDBlib::load >> file "<<filename<<" not found !!!!"<<G4endl;
#else
    cout << "PDBlib::load >> file "<<filename<<" not found !!!!"<<endl;
#endif
  }
  else
  {
    int nbAtomTot=0;        // total number of atoms
    int nbAtom=0;           // total number of atoms in a residue
    int numAtomInRes=0;     // number of an atom in a residue
    int nbResidue=0;        // total number of residues
    int nbMolecule=0;       // total number of molecule

    Atom * AtomicOld = NULL;
    Atom * AtomicNext = NULL;

    int serial;         //Atom serial number
    string atomName;    //Atom name
    string element;     //Element Symbol
    string resName;     //Residue name for this atom
    double x,y,z;        //Orthogonal coordinates in Angstroms
    double occupancy;    //Occupancy

    Residue * residueOld = NULL;
    Residue * residueFirst = NULL;
    Residue * residueNext = NULL;

    Molecule * moleculeFirst = NULL;
    Molecule * moleculeOld = NULL;
    Molecule * moleculeNext = NULL;

    //NEW variable to draw a fitting cylinder if z oriented
    //=> fitting box
    double minGlobZ,maxGlobZ;
    double minGlobX,maxGlobX;
    double minGlobY,maxGlobY;
    double minX,maxX,minY,maxY,minZ,maxZ; //Sort of 'mother volume' box

#ifdef GEANT4
    minGlobZ=-2000000000;
    minGlobX=-2000000000;
    minGlobY=-2000000000;
    maxGlobZ=2000000000;
    maxGlobX=2000000000;
    maxGlobY=2000000000;
    minX=-2000000000;
    minY=-2000000000;
    minZ=-2000000000;
    maxX=2000000000;
    maxY=2000000000;
    maxZ=2000000000;
#else
    minGlobZ=numeric_limits<double>::min();
    minGlobX=numeric_limits<double>::min();
    minGlobY=numeric_limits<double>::min();
    maxGlobZ=numeric_limits<double>::max();
    maxGlobX=numeric_limits<double>::max();
    maxGlobY=numeric_limits<double>::max();
    minX=numeric_limits<double>::min();
    minY=numeric_limits<double>::min();
    minZ=numeric_limits<double>::min();
    maxX=numeric_limits<double>::max();
    maxY=numeric_limits<double>::max();
    maxZ=numeric_limits<double>::max();
#endif

    int lastResSeq=-1;
    int resSeq=0;

    string nameMol;

    unsigned short int numModel=0;
    unsigned short int model=0;

    //Number of TER
    int ter=0;
    //Number of TER (chain) max for this file
#ifdef GEANT4
    int terMax=2000000000;
#else
    int terMax=numeric_limits<int>::max();
#endif
    if (!infile.eof())
    {
      getline(infile, sLine);
      std::size_t found = sLine.find("DNA");
      if (found!=std::string::npos)
      {
        terMax=2;
        isDNA=1;
      }
      else
        isDNA=0;
      //If PDB file have not a header line
      found = sLine.find("HEADER");
      if (found==std::string::npos)
      {
        infile.close();
        infile.open(filename.c_str());
        //Specific to Geant4, use std::cout instead
#ifdef GEANT4
        G4cout<<"PDBlib::load >> No header found !!!!"<<G4endl;
#else
        cout<<"PDBlib::load >> No header found !!!!"<<endl;
#endif
      }

    }

    while (!infile.eof())
    {
      getline(infile, sLine);
      //In the case of several models
      if ((sLine.substr(0,6)).compare("NUMMDL") == 0)
      {
      istringstream ((sLine.substr(10,4))) >> numModel;
      }
      if ((numModel > 0) && ( (sLine.substr(0,6)).compare("MODEL ") == 0))
      {
        istringstream ((sLine.substr(10,4))) >> model;
        if (model > 1) break;
      }
      //For verification of residue sequence
      if ((sLine.substr(0,6)).compare("SEQRES") == 0)
      {
        //Create list of molecule here
#ifdef GEANT4
        if (verbose > 1) G4cout << sLine << G4endl;
#else
        if (verbose > 1) cout << sLine << endl;
#endif
      }

      //Coordinate section
      if ((numModel > 0) && ( (sLine.substr(0,6)).compare("ENDMDL") == 0))
      {;
      }
      else if ((sLine.substr(0,6)).compare("TER   ") == 0) //3 spaces
      {
        //Currently retrieve only the two first chains(TER) => two DNA strands
        ter++;
        if (ter > terMax) break;

        //if (verbose>1) G4cout << sLine << G4endl;
        /************ Begin TER ******************/
        lastResSeq=-1;
        resSeq=0;

        AtomicOld->SetNext(NULL);
        residueOld->SetNext(NULL);
        residueOld->fNbAtom=nbAtom;

        //Molecule creation:
        if (moleculeOld == NULL)
        {
          nameMol=filename;//+numModel
          moleculeOld =new Molecule(nameMol,nbMolecule);
          moleculeOld->SetFirst(residueFirst);
          moleculeOld->fNbResidue=nbResidue;
          moleculeFirst = moleculeOld;
        }
        else
        {
          moleculeNext = new Molecule(nameMol,nbMolecule);
          moleculeOld->SetNext(moleculeNext);
          moleculeNext ->SetFirst(residueFirst);
          moleculeNext ->fNbResidue=nbResidue;
          moleculeOld = moleculeNext;
        }

        nbMolecule++;
        moleculeOld->SetNext(NULL);
        moleculeOld->fCenterX = (int)((minX + maxX) /2.);
        moleculeOld->fCenterY = (int)((minY + maxY) /2.);
        moleculeOld->fCenterZ = (int)((minZ + maxZ) /2.);
        moleculeOld->fMaxGlobZ = maxGlobZ;
        moleculeOld->fMinGlobZ = minGlobZ;
        moleculeOld->fMaxGlobX = maxGlobX;
        moleculeOld->fMinGlobX = minGlobX;
        moleculeOld->fMaxGlobY = maxGlobY;
        moleculeOld->fMinGlobY = minGlobY;

#ifdef GEANT4
        minGlobZ=-2000000000;
        minGlobX=-2000000000;
        minGlobY=-2000000000;
        maxGlobZ=2000000000;
        maxGlobX=2000000000;
        maxGlobY=2000000000;
        minX=-2000000000;
        minY=-2000000000;
        minZ=-2000000000;
        maxX=2000000000;
        maxY=2000000000;
        maxZ=2000000000;
#else
        minGlobZ=numeric_limits<double>::min();
        minGlobX=numeric_limits<double>::min();
        minGlobY=numeric_limits<double>::min();
        maxGlobZ=numeric_limits<double>::max();
        maxGlobX=numeric_limits<double>::max();
        maxGlobY=numeric_limits<double>::max();
        minX=numeric_limits<double>::min();
        minY=numeric_limits<double>::min();
        minZ=numeric_limits<double>::min();
        maxX=numeric_limits<double>::max();
        maxY=numeric_limits<double>::max();
        maxZ=numeric_limits<double>::max();
#endif

        nbAtom=0;
        numAtomInRes=0;
        nbResidue=0;
        AtomicOld = NULL;
        AtomicNext = NULL;
        residueOld = NULL;
        residueFirst = NULL;
        residueNext = NULL;

      }
      else if ((sLine.substr(0,6)).compare("ATOM  ") == 0)
      {

        /************ Begin ATOM ******************/
        //serial
        istringstream ((sLine.substr(6,5))) >> serial;

        //To be improved
        //atomName :
        atomName=sLine.substr(12,4);
        if (atomName.substr(0,1).compare(" ") == 0 )
          element=sLine.substr(13,1);
        else
          element=sLine.substr(12,1);

        // set Van der Waals radius expressed in Angstrom
        double vdwRadius;
        if(element=="H")
        {
          vdwRadius=1.2;
        }
        else if (element=="C")
        {
          vdwRadius=1.7;
        }
        else if (element=="N")
        {
          vdwRadius=1.55;
        }
        else if (element=="O")
        {
          vdwRadius=1.52;
        }
        else if (element=="P")
        {
          vdwRadius=1.8;
        }
        else if (element=="S")
        {
          vdwRadius=1.8;
        }
        else
        {
#ifdef GEANT4
          G4cout << "Element not recognized : " << element << G4endl;
          G4cout << "Stop now" << G4endl;
#else
          cout << "Element not recognized : " << element << endl;
          cout << "Stop now" << endl;
#endif
          exit(1);
        }

        {
          nbAtomTot++;

          //resName :
          resName=sLine.substr(17,3);
          //resSeq :
          istringstream ((sLine.substr(22,4))) >> resSeq;
          //x,y,z :
          stringstream ((sLine.substr(30,8))) >> x;
          stringstream ((sLine.substr(38,8))) >> y;
          stringstream ((sLine.substr(46,8))) >> z;
          //occupancy :
          occupancy=1.;

          if (minGlobZ  < z) minGlobZ=z;
          if (maxGlobZ > z) maxGlobZ=z;
          if (minGlobX  < x) minGlobX=x;
          if (maxGlobX > x) maxGlobX=x;
          if (minGlobY  < y) minGlobY=y;
          if (maxGlobY > y) maxGlobY=y;
          if (minX  > x) minX=x;
          if (maxX < x) maxX=x;
          if (minY  > y) minY=y;
          if (maxY < y) maxY=y;
          if (minZ  > z) minZ=z;
          if (maxZ < z) maxZ=z;

          //treatment for Atoms:
          if (AtomicOld == NULL)
          {
            AtomicOld = new Atom(serial,
                                 atomName,
                                 "",
                                 0,
                                 0,
                                 x,y,z,
                                 vdwRadius,
                                 occupancy,
                                 0,
                                 element);
            AtomicOld->SetNext(NULL);//If only one Atom inside residue
          }
          else
          {
            AtomicNext = new Atom(serial,
                                  atomName,
                                  "",
                                  0,
                                  0,
                                  x,y,z,
                                  vdwRadius,
                                  occupancy,
                                  0,
                                  element);
            AtomicOld->SetNext(AtomicNext);
            AtomicOld = AtomicNext;
          }
          nbAtom++;
        }//END if (element!="H")
        /****************************Begin Residue************************/
        //treatment for residues:
        if (residueOld == NULL)
        {
#ifdef GEANT4
          if (verbose>2) G4cout << "residueOld == NULL"<<G4endl;
#else
          if (verbose>2) cout << "residueOld == NULL"<<endl;
#endif
          AtomicOld->fNumInRes=0;
          residueOld = new Residue(resName,resSeq);
          residueOld->SetFirst(AtomicOld);
          residueOld->SetNext(NULL);
          residueFirst = residueOld;
          lastResSeq=resSeq;
          nbResidue++;
        }
        else
        {
          if (lastResSeq==resSeq)
          {
            numAtomInRes++;
            AtomicOld->fNumInRes=numAtomInRes;
          }
          else
          {
            numAtomInRes=0;
            AtomicOld->fNumInRes=numAtomInRes;

            residueNext = new Residue(resName,resSeq);
            residueNext->SetFirst(AtomicOld);
            residueOld->SetNext(residueNext);
            residueOld->fNbAtom=nbAtom-1;

            nbAtom=1;
            residueOld = residueNext;
            lastResSeq=resSeq;
            nbResidue++;
          }
        }
      }//END if Atom
    }//END while (!infile.eof())

    infile.close();
    return moleculeFirst;
  }//END else if (!infile)
  return NULL;
}

Here is the call graph for this function:

Here is the caller graph for this function:

Member Data Documentation

fNbNucleotidsPerStrand

int PDBlib::fNbNucleotidsPerStrand

private

Number of nucleotid per strand.

Definition at line 97 of file PDBlib.hh.

---

The documentation for this class was generated from the following files:

• PDBlib.hh
• PDBlib.cc