replace the PDB res nums with atom positions:

@param afpChain
@param ca1
@param ca2
17926
public static void rebuildAFPChain(AFPChain afpChain, Atom[] ca1, Atom[] ca2){

		if ( afpChain.getAlgorithmName() == null) {
			afpChain.setAlgorithmName(DEFAULT_ALGORITHM_NAME);
		}
		if ( afpChain.getVersion() == null){
			afpChain.setVersion("1.0");
		}

		int blockNum  = afpChain.getBlockNum();
		int ca1Length = afpChain.getCa1Length();
		int ca2Length = afpChain.getCa2Length();

		int minLength = Math.min(ca1Length, ca2Length);
		int[][][] optAln = new int[blockNum][2][minLength];

		int[][][] blockResList = afpChain.getBlockResList();
		if ( blockResList == null){
			blockResList = new int[blockNum][2][minLength];
		}
		int[] optLen = afpChain.getOptLen();

		String[][][] pdbAln = afpChain.getPdbAln();
		int[] verifiedOptLen = null;
		if ( optLen != null)
		  verifiedOptLen = afpChain.getOptLen().clone();
		else {
			logger.warn("did not find optimal alignment, building up empty alignment.");
			optLen = new int[1];
			optLen[0] = 0;
		}
		for (int blockNr = 0 ; blockNr < blockNum ; blockNr++){

			//System.out.println("got block " + blockNr + " size: " + optLen[blockNr]);
			int verifiedEQR = -1;
			for ( int eqrNr = 0 ; eqrNr < optLen[blockNr] ; eqrNr++ ){
				String pdbResnum1 = pdbAln[blockNr][0][eqrNr];
				String pdbResnum2 = pdbAln[blockNr][1][eqrNr];

				//System.out.println(blockNr + " " + eqrNr + " got resnum: " + pdbResnum1 + " " + pdbResnum2);
				String[] spl1 = pdbResnum1.split(":");
				String[] spl2 = pdbResnum2.split(":");

				String chain1 = spl1[0];
				String pdbres1 = spl1[1];

				String chain2 = spl2[0];
				String pdbres2 = spl2[1];

				int pos1 = getPositionForPDBresunm(pdbres1,chain1,ca1);
				int pos2 = getPositionForPDBresunm(pdbres2,chain2,ca2);

				if ( pos1 == -1 || pos2 == -1 ){
					// this can happen when parsing old files that contained Calcium atoms...
					logger.warn("pos1: {} (residue {}), pos2: {} (residue {}), should never be -1. Probably parsing an old file.",
							pos1, pdbResnum1, pos2, pdbResnum2);
					verifiedOptLen[blockNr]-- ;
					continue;
				}

				verifiedEQR++;
				//System.out.println(blockNr + " " + eqrNr + " " + pos1 + " " + pos2);
				optAln[blockNr][0][verifiedEQR] = pos1;
				optAln[blockNr][1][verifiedEQR] = pos2;
				blockResList[blockNr][0][verifiedEQR] = pos1;
				blockResList[blockNr][1][verifiedEQR] = pos2;
			}
		}

		afpChain.setOptLen(verifiedOptLen);
		afpChain.setOptAln(optAln);
		afpChain.setBlockResList(blockResList);
		// build up alignment image:
		AFPAlignmentDisplay.getAlign(afpChain, ca1, ca2);


	}
1060642
private static void permuteAFPChain(AFPChain afpChain, int cp) {

		int ca2len = afpChain.getCa2Length();

		//fix up cp to be positive
		if(cp == 0) {
			return;
		}
		if(cp < 0) {
			cp = ca2len+cp;
		}
		if(cp < 0 || cp >= ca2len) {
			throw new ArrayIndexOutOfBoundsException(
					"Permutation point ("+cp+") must be between -ca2.length and ca2.length-1" );
		}

		// Fix up optAln
		permuteOptAln(afpChain,cp);

		if(afpChain.getBlockNum() > 1)
			afpChain.setSequentialAlignment(false);
		// fix up matrices
		// ca1 corresponds to row indices, while ca2 corresponds to column indices.
		afpChain.setDistanceMatrix(permuteMatrix(afpChain.getDistanceMatrix(),0,-cp));
		// this is square, so permute both
		afpChain.setDisTable2(permuteMatrix(afpChain.getDisTable2(),-cp,-cp));

		//TODO fix up other AFP parameters?

	}
1178349
private static int  calcGap(AFP afp1 , AFP afp2)
	{
		int     g = (afp1.getP1() - afp2.getP1()) - (afp1.getP2() - afp2.getP2());
		if(g < 0)       g = -g;
		return g;
	}
446107
public static void optimizeAln(FatCatParameters params, AFPChain afpChain,Atom[] ca1, Atom[] ca2) throws StructureException
	{

		int minLen = afpChain.getMinLen();
		int fragLen = params.getFragLen();


		long optStart = System.currentTimeMillis();
		int     i, a, k, p1, p2, bk, b1, b2, e1, e2, a1, a2;
		int     iniLen;
		int[][]     iniSet = new int[2][minLen];
		int     maxi = 100;

		int[][][] optAln = afpChain.getOptAln();
		int[] optLen = afpChain.getOptLen();

		int maxTra = params.getMaxTra();
		double[] optRmsd = afpChain.getOptRmsd();

		int blockNum = afpChain.getBlockNum();

		if(optAln == null)      {
			optAln     = new int[maxTra+1][2][minLen];
			optLen     = new int[maxTra+1];
			afpChain.setOptLen(optLen);
			optRmsd    = new double[maxTra+1];
			afpChain.setOptRmsd(optRmsd);
		}

		List<AFP> afpSet = afpChain.getAfpSet();

		int optLength         = afpChain.getOptLength();
		int[] afpChainList    = afpChain.getAfpChainList();
		int[] block2Afp       = afpChain.getBlock2Afp();
		int[] blockSize       = afpChain.getBlockSize();


		if (debug)
			System.out.println("AFPOptimizer got blockNum: " +blockNum);
		//optimize each alignment defined by a block
		b1 = b2 = e1 = e2 = optLength = 0;
		for(bk = 0; bk < blockNum; bk ++)       {
			//initial aligned position
			iniLen = 0;
			if(bk > 0)      {
				b1 = e1;
				b2 = e2;
			}
			if(bk < blockNum - 1)   {
				a1 = afpChainList[block2Afp[bk] + blockSize[bk] - 1]; //the last AFP in current block
				a2 = afpChainList[block2Afp[bk + 1]];  //the first AFP in next block
				e1 = (afpSet.get(a1).getP1() + fragLen +  afpSet.get(a2).getP1()) / 2;
				e2 = (afpSet.get(a1).getP2() + fragLen +  afpSet.get(a2).getP2()) / 2;
			} //use the middle point of the current and next AFPs. old (starting point of next AFP)
			else    {
				e1 = ca1.length;
				e2 = ca2.length;
			}


			for(i = block2Afp[bk]; i < block2Afp[bk] + blockSize[bk]; i ++) {
				a = afpChainList[i];
				p1 = afpSet.get(a).getP1();
				p2 = afpSet.get(a).getP2();
				for(k = 0; k < afpSet.get(a).getFragLen(); k ++)     {
					iniSet[0][iniLen] = p1 + k - b1; //note -b1
					iniSet[1][iniLen] = p2 + k - b2; //note -b2
					iniLen ++;

				}
			}
			//optimize the align by dynamic programming & constraint the optimization region
			if(debug) {
				System.err.println(String.format("optimize block %d (%d afp), region %d-%d(len %d), %d-%d(len %d)\n",
						bk, blockSize[bk], b1, e1, e1-b1, b2, e2, e2-b2));

				System.err.println(" initial alignment Length: " + iniLen );
			}

			StructureAlignmentOptimizer opt = new StructureAlignmentOptimizer(b1,e1, ca1, b2,e2, ca2, iniLen, iniSet);
			opt.runOptimization(maxi);
			optRmsd[bk] = opt.optimizeResult(optLen,bk,optAln[bk]);

			//System.out.println(optRmsd[bk]);
			// SALNOPT *opt = new SALNOPT(e1-b1, &pro1->caCod[3 * b1], e2-b2, &pro2->caCod[3 * b2], iniLen, iniSet, maxi);
			// optRmsd[bk] = opt->OptimizeResult(&optLen[bk], optAln[bk]);

			if(debug)
				System.out.println(String.format(" optimized len=%d, rmsd %f\n", optLen[bk], optRmsd[bk]));

			for(i = 0; i < optLen[bk]; i ++)        {
				optAln[bk][0][i] += b1; //restore the position
				optAln[bk][1][i] += b2; //restore the position
			}
			optLength += optLen[bk];


		}

		long optEnd = System.currentTimeMillis();
		if(debug)       System.out.println("complete AlignOpt " + (optEnd-optStart) +"\n");

		if(optLength < minLen) {
			int[][][] optAln_trim = new int[maxTra + 1][2][optLength];
			for (i = 0; i < maxTra + 1; i ++) {
				System.arraycopy(optAln[i][0], 0, optAln_trim[i][0], 0, optLength);
				System.arraycopy(optAln[i][1], 0, optAln_trim[i][1], 0, optLength);
			}
			afpChain.setOptAln(optAln_trim);
		} else {
			afpChain.setOptAln(optAln);
		}

		afpChain.setBlockNum(blockNum);
		afpChain.setOptLength(optLength);
		afpChain.setAfpChainList(afpChainList);
		afpChain.setBlock2Afp(block2Afp);
		afpChain.setBlockSize(blockSize);


	}
332217
private static int getPositionForPDBresunm(String pdbresnum, String authId , Atom[] atoms){
		ResidueNumber residueNumber =  ResidueNumber.fromString(pdbresnum);
		residueNumber.setChainName(authId);

		boolean blankChain = authId == null || authId.equalsIgnoreCase("null") || authId.equals("_");

		for ( int i =0; i< atoms.length ;i++){
			Group g = atoms[i].getGroup();

			// match _ to any chain
			if( blankChain ) {
				residueNumber.setChainName(g.getChain().getName());
			}

			//System.out.println(g.getResidueNumber() + "< ? >" + residueNumber +"<");
			if ( g.getResidueNumber().equals(residueNumber)){
				//System.out.println(g + " == " + residueNumber );
				Chain c = g.getChain();
				if ( blankChain || c.getName().equals(authId)){
					return i;
				}
			}
		}
		return -1;
	}
142830
public void align(Structure s1, String chainId1, Structure s2,
			String chainId2, StrucAligParameters params)
			throws StructureException {
		reset();
		this.params = params;

		Chain c1 = s1.getPolyChainByPDB(chainId1);
		Chain c2 = s2.getPolyChainByPDB(chainId2);

		Structure s3 = new StructureImpl();
		s3.addChain(c1);

		Structure s4 = new StructureImpl();
		s4.addChain(c2);

		Atom[] ca1 = getAlignmentAtoms(s3);
		Atom[] ca2 = getAlignmentAtoms(s4);

		notifyStartingAlignment(s1.getName(), ca1, s2.getName(), ca2);
		align(ca1, ca2, params);
	}
333700
private void calculateHBonds() {
		/**
		 * More efficient method for calculating C-Alpha pairs
		 */
		if (groups.length < 5) return;
		Iterator<AtomContact> otu = contactSet.iterator();
		while(otu.hasNext()){
			AtomContact ac = otu.next();
			Pair<Atom> pair = ac.getPair();
			Group g1 = pair.getFirst().getGroup();
			Group g2 = pair.getSecond().getGroup();
			// Now I need to get the index of the Group in the list groups
			int i = indResMap.get(g1.getResidueNumber());
			int j = indResMap.get(g2.getResidueNumber());
			// Now check this
			checkAddHBond(i,j);
			//"backwards" hbonds are not allowed
			if (j!=(i+1)) checkAddHBond(j,i);
		}
	}
125799
def match_RCSB_pdb_chains(pdb_id1, pdb_id2, cut_off = 60.0, allow_multiple_matches = False, multiple_match_error_margin = 3.0, use_seqres_sequences_if_possible = True, strict = True):
    ''''''
    try:
        stage = pdb_id1
        pdb_1 = PDB(retrieve_pdb(pdb_id1), strict = strict)
        stage = pdb_id2
        pdb_2 = PDB(retrieve_pdb(pdb_id2), strict = strict)
    except (PDBParsingException, NonCanonicalResidueException, PDBValidationException), e:
        raise PDBParsingException("An error occurred while loading %s: '%s'" % (stage, str(e)))

    return match_pdb_chains(pdb_1, pdb_id1, pdb_2, pdb_id2, cut_off = cut_off, allow_multiple_matches = allow_multiple_matches, multiple_match_error_margin = multiple_match_error_margin, use_seqres_sequences_if_possible = use_seqres_sequences_if_possible)
434350
public static void updateScore(FatCatParameters params, AFPChain afpChain)
	{
		int     i, j, bknow, bkold, g1, g2;


		afpChain.setConn(0d);
		afpChain.setDVar(0d);

		int blockNum = afpChain.getBlockNum();
		int alignScoreUpdate = 0;
		double[] blockScore = afpChain.getBlockScore();
		int[] blockGap = afpChain.getBlockGap();
		int[] blockSize =afpChain.getBlockSize();
		int[] afpChainList = afpChain.getAfpChainList();
		List<AFP>afpSet = afpChain.getAfpSet();
		int[] block2Afp = afpChain.getBlock2Afp();

		double torsionPenalty = params.getTorsionPenalty();


		bkold = 0;
		for(i = 0; i < blockNum; i ++)  {
			blockScore[i] = 0;
			blockGap[i] = 0;
			for(j = 0; j < blockSize[i]; j ++)      {
				bknow = afpChainList[block2Afp[i] + j];
				if(j == 0)      {
					blockScore[i] = afpSet.get(bknow).getScore();
				}
				else    {
					AFPChainer.afpPairConn(bkold, bknow, params, afpChain); //note: j, i
					Double conn = afpChain.getConn();
					blockScore[i] += afpSet.get(bknow).getScore() + conn;
					g1 = afpSet.get(bknow).getP1() - afpSet.get(bkold).getP1() - afpSet.get(bkold).getFragLen();
					g2 = afpSet.get(bknow).getP2() - afpSet.get(bkold).getP2() - afpSet.get(bkold).getFragLen();
					blockGap[i] += (g1 > g2)?g1:g2;
				}
				bkold = bknow;
			}
			alignScoreUpdate += blockScore[i];
		}
		if(blockNum >= 2)       {
			alignScoreUpdate += (blockNum - 1) * torsionPenalty;
		}

		afpChain.setBlockGap(blockGap);
		afpChain.setAlignScoreUpdate(alignScoreUpdate);
		afpChain.setBlockScore(blockScore);
		afpChain.setBlockSize(blockSize);
		afpChain.setAfpChainList(afpChainList);
		afpChain.setBlock2Afp(block2Afp);
	}
40098
public static void shiftCA2(AFPChain afpChain, Atom[] ca2,  Matrix m, Atom shift, Group[] twistedGroups) {

		int i = -1;
		for (Atom a: ca2){
			i++;
			Group g = a.getGroup();

			Calc.rotate(g,m);
			Calc.shift(g, shift);

			if (g.hasAltLoc()){
				for (Group alt: g.getAltLocs()){
					for (Atom alta : alt.getAtoms()){
						if ( g.getAtoms().contains(alta))
							continue;
						Calc.rotate(alta,m);
						Calc.shift(alta,shift);
					}
				}
			}
			twistedGroups[i]=g;
		}
	}
1175569
def get_differing_atom_residue_ids(self, pdb_name, pdb_list):
        ''''''

        assert(pdb_name in self.pdb_names)
        assert(set(pdb_list).intersection(set(self.pdb_names)) == set(pdb_list)) # the names in pdb_list must be in pdb_names

        differing_atom_residue_ids = set()
        for other_pdb in pdb_list:
            differing_atom_residue_ids = differing_atom_residue_ids.union(set(self.differing_atom_residue_ids[(pdb_name, other_pdb)]))

        return sorted(differing_atom_residue_ids)
897270
private static void putToMapIfNotThere(Map<String, Integer> symbolsMap,
                                           String text, int sid)
    {
        // When there's a duplicate mapping for the symbol text, don't
        // replace the lower sid. This pattern avoids double-lookup in the
        // common scenario and only performs a second lookup when there's
        // a duplicate.
        Integer extantSid = symbolsMap.put(text, sid);
        if (extantSid != null)
        {
            // We always insert symbols with increasing sids
            assert extantSid < sid;
            symbolsMap.put(text, extantSid);
        }
    }
678223
public static AFPChain fastaToAfpChain(SequencePair<Sequence<AminoAcidCompound>, AminoAcidCompound> alignment,
			Structure structure1, Structure structure2) throws StructureException {
		List<AlignedSequence<Sequence<AminoAcidCompound>, AminoAcidCompound>> seqs = alignment.getAlignedSequences();
		StringBuilder sb1 = new StringBuilder();
		for (AminoAcidCompound a : seqs.get(0)) {
			sb1.append(a.getBase());
		}
		try {
			ProteinSequence seq1 = new ProteinSequence(sb1.toString());
			StringBuilder sb2 = new StringBuilder();
			for (AminoAcidCompound a : seqs.get(1)) {
				sb1.append(a.getBase());
			}
			ProteinSequence seq2 = new ProteinSequence(sb2.toString());
			LinkedHashMap<String, ProteinSequence> map = new LinkedHashMap<String, ProteinSequence>();
			map.put(structure1.getName(), seq1);
			map.put(structure2.getName(), seq2);
			return fastaToAfpChain(map, structure1, structure2);
		} catch (CompoundNotFoundException e) {
			logger.error("Unexpected error while creating protein sequences: {}. This is most likely a bug.",e.getMessage());
			return null;
		}
	}
943620
public static int getBlockNrForAlignPos(AFPChain afpChain, int aligPos){
		// moved here from DisplayAFP;

		int blockNum = afpChain.getBlockNum();

		int[] optLen = afpChain.getOptLen();
		int[][][] optAln = afpChain.getOptAln();

		int len = 0;
		int p1b=0;
		int p2b=0;

		for(int i = 0; i < blockNum; i ++)  {

			for(int j = 0; j < optLen[i]; j ++) {

				int p1 = optAln[i][0][j];
				int p2 = optAln[i][1][j];

				if (len != 0) {
					// check for gapped region
					int lmax = (p1 - p1b - 1)>(p2 - p2b - 1)?(p1 - p1b - 1):(p2 - p2b - 1);
					for(int k = 0; k < lmax; k ++)      {
						len++;
					}
				}


				p1b = p1;
				p2b = p2;
				if ( len >= aligPos) {

					return i;
				}
				len++;
			}
		}

		return blockNum;

	}
200715
private void rotateShiftAtoms(Atom[] ca){

		for (int i  = 0 ; i < ca.length; i++){
			Atom c = ca[i];

			Calc.rotate(c,currentRotMatrix);
			Calc.shift(c,currentTranMatrix);
			//System.out.println("after " + c);
			ca[i] = c;
		}
		//System.out.println("after " + ca[0]);
	}
1069575
def update_atom_numbers(self):
        
        self._ctab_data['first_atom_number'] = self.first_atom.atom_number
        self._ctab_data['second_atom_number'] = self.second_atom.atom_number
173681
public AFPChain alignPermuted(Atom[] ca1, Atom[] ca2, Object param, int cp) throws StructureException {
		// initial permutation
		permuteArray(ca2,cp);

		// perform alignment
		AFPChain afpChain = super.align(ca1, ca2, param);

		// un-permute alignment
		permuteAFPChain(afpChain, -cp);

		if(afpChain.getName2() != null) {
			afpChain.setName2(afpChain.getName2()+" CP="+cp);
		}

		// Specify the permuted
		return afpChain;
	}