--- trunk/OOPSE/libBASS/MoleculeStamp.cpp	2004/06/01 21:44:54	1216
+++ trunk/OOPSE/libBASS/MoleculeStamp.cpp	2004/06/10 14:59:14	1256
@@ -308,6 +308,8 @@ char* MoleculeStamp::addRigidBody( RigidBodyStamp* the
 char* MoleculeStamp::addRigidBody( RigidBodyStamp* the_rigidbody, 
                                    int rigidBodyIndex ){
   
+
+  printf("rigidBodyIndex = %d\n", rigidBodyIndex);
   if( have_rigidbodies && rigidBodyIndex < n_rigidbodies ) 
     rigidBodies[rigidBodyIndex] = the_rigidbody;
   else {
@@ -454,3 +456,81 @@ char* MoleculeStamp::checkMe( void ){
 
   return NULL;
 }  
+
+
+//Function Name: isBondInSameRigidBody
+//Return true is both atoms of the bond belong to the same rigid body, otherwise return false
+bool MoleculeStamp::isBondInSameRigidBody(BondStamp* bond){
+  int rbA;
+  int rbB;
+  int consAtomA;
+  int consAtomB;
+
+  return isAtomInRigidBody(bond->getA(),rbA, consAtomA) &&
+                isAtomInRigidBody(bond->getB(),rbB, consAtomB);
+}
+
+
+// Function Name isAtomInRigidBody 
+//return false if atom does not belong to a rigid body otherwise return true and set whichRigidBody 
+//and consAtomIndex
+//atomIndex : the index of atom in component
+//whichRigidBody: the index of rigidbody in component
+//consAtomIndex:  the position of joint atom apears in  rigidbody's definition
+bool MoleculeStamp::isAtomInRigidBody(int atomIndex, int& whichRigidBody, int& consAtomIndex){
+  RigidBodyStamp* rbStamp;
+  int numRb;
+  int numAtom;
+
+  numRb = this->getNRigidBodies();
+  
+  for(int i = 0 ; i < numRb; i++){
+    rbStamp = this->getRigidBody(i);
+    numAtom = rbStamp->getNMembers();
+    for(int j = 0; j < numAtom; j++)
+      if (rbStamp->getMember(j) == atomIndex){
+        whichRigidBody = i;
+        consAtomIndex = j;
+        return true;
+      }
+  }
+
+  return false;
+   
+}
+
+//return the position of joint atom apears in  rigidbody's definition
+//for the time being, we will use the most inefficient algorithm, the complexity is O(N2)
+//actually we could improve the complexity to O(NlgN) by sorting the atom index in rigid body first
+vector<pair<int, int> > MoleculeStamp::getJointAtoms(int rb1, int rb2){
+  RigidBodyStamp* rbStamp1;
+  RigidBodyStamp* rbStamp2;
+  int natomInRb1;
+  int natomInRb2;
+  int atomIndex1;
+  int atomIndex2;
+  vector<pair<int, int> > jointAtomIndexPair;
+  
+  rbStamp1 = this->getRigidBody(rb1);
+  natomInRb1 =rbStamp1->getNMembers();
+
+  rbStamp2 = this->getRigidBody(rb2);
+  natomInRb2 =rbStamp2->getNMembers();
+
+  for(int i = 0; i < natomInRb1; i++){
+    atomIndex1 = rbStamp1->getMember(i);
+      
+    for(int j= 0; j < natomInRb1; j++){
+      atomIndex2 = rbStamp2->getMember(j);
+
+      if(atomIndex1 == atomIndex2){
+        jointAtomIndexPair.push_back(make_pair(i, j));
+        break;
+      }
+      
+    }//end for(j =0)
+
+  }//end for (i = 0)
+
+  return jointAtomIndexPair;
+}