#include "BiopolymerClass.h"
 #ifndef ResidueStretchContainer_H_
 #define ResidueStretchContainer_H_

class ResidueStretchContainer{
	private:
	vector<ResidueStretch> residueStretchVector;
	public:
	ResidueStretchContainer() {};

	void clear(){residueStretchVector.clear();};
        vector<ResidueStretch> getResidueStretchVector(){return residueStretchVector;};
	vector <MobilizerStretch> getMobilizerStretchVector(String bondMobilityString) {
		vector <MobilizerStretch> myMobilizerStretchVector;
		for (int i = 0; i < (int)residueStretchVector.size(); i++){
                        MobilizerStretch myMobilizerStretch = residueStretchVector[i].mobilizerStretch(bondMobilityString);
			myMobilizerStretchVector.push_back(myMobilizerStretch);// residueStretchVector[i].mobilizerStretch(bondMobilityString));
              
                        cout<<__FILE__<<":"<<__LINE__<<" For chain "<<myMobilizerStretch.Chain<< " start res "<<myMobilizerStretch.StartResidue.outString()<<" end res "<< myMobilizerStretch.EndResidue.outString() <<" mobility " << myMobilizerStretch.BondMobility<<endl;//" and chain "<< targetChain<<" residue "<<targetResidue.outString()<<" distance is "<<myDistance<<" Å "<<endl;
		}
		return myMobilizerStretchVector;
	};
	void addResidueStretchToVector(ResidueStretch residueStretch){residueStretchVector.push_back(residueStretch);};
	bool residuesAreWithinCutoff(double radius, BiopolymerClass & referenceBiopolymerClass, ResidueID referenceResidue, BiopolymerClass & targetBiopolymerClass, ResidueID targetResidue, BiopolymerClassContainer myBiopolymerClassContainer){
                String referenceChain = referenceBiopolymerClass.getChainID();
                String targetChain = targetBiopolymerClass.getChainID();
		double myDistance = residueDistance(referenceBiopolymerClass, referenceResidue, targetBiopolymerClass, targetResidue, myBiopolymerClassContainer);
		if (myDistance <= radius) {return true;} else {return false;}
		
	};
				

	double  residueDistance( BiopolymerClass & referenceBiopolymerClass, ResidueID referenceResidue, BiopolymerClass & targetBiopolymerClass, ResidueID targetResidue, BiopolymerClassContainer myBiopolymerClassContainer){
                String referenceChain = referenceBiopolymerClass.getChainID();
                String targetChain = targetBiopolymerClass.getChainID();
		double myDistance = (double)(
	                    	targetBiopolymerClass.calcDefaultAtomLocationInGroundFrame  ( targetResidue, targetBiopolymerClass.getRepresentativeAtomName())
                    		- referenceBiopolymerClass.calcDefaultAtomLocationInGroundFrame(referenceResidue, referenceBiopolymerClass.getRepresentativeAtomName())
                    	).norm()*10; // convert to Å
                cout<<__FILE__<<":"<<__LINE__<<" For chain "<<referenceChain<<" residue "<<referenceResidue.outString()<<" and chain "<< targetChain<<" residue "<<targetResidue.outString()<<" distance is "<<myDistance<<" Å "<<endl;
                return myDistance;
		
	};
				




	void addTargetChainResidue(double radius, BiopolymerClass & referenceBiopolymerClass, ResidueID referenceResidue, BiopolymerClass & targetBiopolymerClass, ResidueID targetResidue,BiopolymerClassContainer myBiopolymerClassContainer){
                String referenceChain = referenceBiopolymerClass.getChainID(); 
                String targetChain = targetBiopolymerClass.getChainID(); 
		if (referenceBiopolymerClass.getChainID().compare(targetBiopolymerClass.getChainID()) == 0) {
			cout<<__FILE__<<":"<<__LINE__<<" Chain "<<referenceChain<<" can't have an interface with itself!"<<endl;
			exit(1);
		};
			residueStretchVector.push_back(ResidueStretch(targetChain, targetResidue, targetResidue));
                
		cout<<__FILE__<<":"<<__LINE__<<" Chain "<<targetChain<<" residue "<<targetResidue.outString()<<" added to interface mobility zone. "<<endl;
	};


	void addMutualChainResidues(double radius, BiopolymerClass & referenceBiopolymerClass, BiopolymerClass & targetBiopolymerClass,BiopolymerClassContainer myBiopolymerClassContainer){
		//BiopolymerClass targetBiopolymerClass = myBiopolymerClassContainer.updBiopolymerClass(targetChain);
		//BiopolymerClass referenceBiopolymerClass = myBiopolymerClassContainer.updBiopolymerClass(referenceChain);

		for (ResidueID referenceResidueID = referenceBiopolymerClass.getFirstResidueID(); referenceResidueID <= referenceBiopolymerClass.getLastResidueID(); referenceBiopolymerClass.incrementResidueID(referenceResidueID)) {
			for (ResidueID targetResidueID = targetBiopolymerClass.getFirstResidueID(); targetResidueID <= targetBiopolymerClass.getLastResidueID(); targetBiopolymerClass.incrementResidueID(targetResidueID))
			{
				addTargetChainResidue(radius, referenceBiopolymerClass, referenceResidueID, targetBiopolymerClass,        targetResidueID,       myBiopolymerClassContainer );	
                                addTargetChainResidue(radius, targetBiopolymerClass,        targetResidueID,       referenceBiopolymerClass, referenceResidueID, myBiopolymerClassContainer );
                                if (targetResidueID == targetBiopolymerClass.getLastResidueID()) {break;  }
			}
                        if (referenceResidueID == referenceBiopolymerClass.getLastResidueID()) {break;  }
                }
	};

	void addAllMutualChainResidues(double radius, String referenceChain, BiopolymerClassContainer myBiopolymerClassContainer){
		BiopolymerClass referenceBiopolymerClass = myBiopolymerClassContainer.updBiopolymerClass(referenceChain);

                cout<<__FILE__<<":"<<__LINE__<<" Interface cutoff  is : "<<radius<<endl;
		for (int i = 0; i < myBiopolymerClassContainer.getNumBiopolymers(); i++)
			if (referenceChain.compare(myBiopolymerClassContainer.updBiopolymerClass(i).getChainID()) != 0) // don't check referenceChain against itself!
			{
 				BiopolymerClass targetBiopolymerClass = myBiopolymerClassContainer.updBiopolymerClass(i);
                                String targetChain = targetBiopolymerClass.getChainID();
				for (ResidueID referenceResidueID = referenceBiopolymerClass.getFirstResidueID(); referenceResidueID <= referenceBiopolymerClass.getLastResidueID(); referenceBiopolymerClass.incrementResidueID(referenceResidueID)){
                                        double lowestDistanceForReferenceResidueID = 11111;
					for (ResidueID targetResidueID = targetBiopolymerClass.getFirstResidueID(); targetResidueID <= targetBiopolymerClass.getLastResidueID(); targetBiopolymerClass.incrementResidueID(targetResidueID)) {
                                                double myDistance = residueDistance(referenceBiopolymerClass,  referenceResidueID, targetBiopolymerClass,  targetResidueID, myBiopolymerClassContainer);
                                                if (myDistance < lowestDistanceForReferenceResidueID) {lowestDistanceForReferenceResidueID = myDistance; }
        					if (myDistance <= radius)
						{
							addTargetChainResidue(radius,referenceBiopolymerClass,referenceResidueID,targetBiopolymerClass,   targetResidueID,   myBiopolymerClassContainer );
							addTargetChainResidue(radius,targetBiopolymerClass,   targetResidueID,   referenceBiopolymerClass,referenceResidueID,myBiopolymerClassContainer );
						} else if (myDistance >  (radius + (3.8 *10))) { // if we can afford to advance 5 residues instead of just 1 ..
							for (int j = 0; j < 9; j++) { // go ahead and skip 9
							    if (targetResidueID == targetBiopolymerClass.getLastResidueID()) {break;  }
							    targetBiopolymerClass.incrementResidueID(targetResidueID); // skip one residue
                                                        }
                                                        
						} else if (myDistance >  (radius + (3.8 * 5))) { // if we can afford to advance 5 residues instead of just 1 ..
							for (int j = 0; j < 4; j++) { // go ahead and skip 4
							    if (targetResidueID == targetBiopolymerClass.getLastResidueID()) {break;  }
							    targetBiopolymerClass.incrementResidueID(targetResidueID); // skip one residue
                                                        }
                                                        
       						}
						else if (myDistance >  (radius + (3.8 * 2))) { // 3.8Å is the per-residue advance for an extended chain.if we can afford to advance 2 residues instead of just 1 ..
                                                        cout<<__FILE__<<":"<<__LINE__<<" Distance is "<<myDistance<<", we can skip a residue"<<endl;
							for (int j = 0; j < 1; j++) { // go ahead and skip 1
							    if (targetResidueID == targetBiopolymerClass.getLastResidueID()) {break;  }
							    targetBiopolymerClass.incrementResidueID(targetResidueID); // skip one residue
                                                        }
                                                }        
						if (targetResidueID == targetBiopolymerClass.getLastResidueID()) {break;  }
                   			}

					int 	lowestPossibleSequenceDistance = int((lowestDistanceForReferenceResidueID-radius)/3.8);
					if (lowestPossibleSequenceDistance >0) {
                                                cout<<__FILE__<<":"<<__LINE__<<" Minimum distance for reference chain is "<<lowestDistanceForReferenceResidueID <<", we can skip some residues"<<endl;
						for (int k = 0 ; k < (lowestPossibleSequenceDistance-2); k++){
							if (referenceResidueID == referenceBiopolymerClass.getLastResidueID()) {break;  }
							referenceBiopolymerClass.incrementResidueID(referenceResidueID); // skip one residue
						}
					}
                                        if (referenceResidueID == referenceBiopolymerClass.getLastResidueID()) {break;  }
                                }
			}
		};
	
};

#endif