public void setupNHVectors() {
    	// First create a map from residue number to the index of corresponding HN atom
    	HashMap<Integer, Integer> HAtoms = new HashMap<Integer, Integer>();
    	PDBData atom;
    	for(int i=0; i<atoms.size(); i++) {
    		atom = atoms.get(i);
    		if(atom.atomName.trim().equals("HN"))
    			HAtoms.put(atom.resID, i);
    	}
    	
    	NHVectors = new ArrayList<NHVector>();
    	Vector bond, n1, n2;
    	for(int i=0; i<atoms.size(); i++) {
    		atom = atoms.get(i);
    		// If we find N atom and there is a HN atom on the same residue then we save the
    		// corresponding NH vector in the list
    		if(atom.atomName.trim().equals("N") && HAtoms.containsKey(atom.resID)) {
    			int Hindex = HAtoms.get(atom.resID);
    			PDBData Hatom = atoms.get(Hindex);
    			NHVector vec = new NHVector();
    			//Save the coordinates of the NH vector
    			bond = new Vector(Hatom.x - atom.x, Hatom.y - atom.y, Hatom.z - atom.z);
    			vec.bond = bond;
    			
     			// Save the coordinates of the N atom
    			vec.Ncoords = new Vector(atom.x, atom.y, atom.z);
    			
     			// Save the indexes of the N and HN atoms
    			vec.indexH = Hindex;
    			vec.indexN = i;
    			vec.resID = atom.resID;
    			
    			// Set the vectors perpendicular to the bond vector. 
    			// Make their lengths equal to that of the bond vector, to simplify future calculations
    			float length = bond.norm(); // bond length
    			float x = bond.x, y = bond.y, z = bond.z;  // Convenience variables to shorten formulas below
    			
    			if(y == 0.0f && z == 0.0f) {  // First deal with the special case to avoid devision by zero
    				n1 = new Vector(0, length, 0);
    				n2 = new Vector(0, 0, length);
    			} else {
        			n1 = new Vector(0, z, -y);
        			n1.multiply(length / n1.norm()); // Set the length of n1 properly
        			
        			n2 = new Vector(-y*y - z*z, x*y, x*z);
        			n2.multiply(length / n2.norm()); // Set the length of n2 properly
    			}
    			
    		//	System.out.printf("Norms: bond=%.3f, n1=%.3f, n2=%.3f\n", bond.norm(), n1.norm(), n2.norm());
    		//	System.out.printf("Dot products: n1*v=%.3f, n2*v=%.3f, n1*n2=%.3f\n", n1.dot(bond), n2.dot(bond), n1.dot(n2));
    			vec.norm1 = n1;
    			vec.norm2 = n2;
    			
    			NHVectors.add(vec);
    		}
    	}
    	
    	System.out.println(NHVectors.size() + " NH vectors were found");
     }