Returns a mapping chain_of_pdb_1 -> List[PDBMutationPair] representing the mutations needed to transform each chain of pdb_1 into the respective chains of pdb_2.
           This function also sets self.seqres_sequence_maps, a mapping Tuple(pdb_id_1, chain ID in pdb_id_1) -> Tuple(pdb_id_2, chain ID in pdb_id_2) -> a SequenceMap representing the mapping of residues between the chains based on the alignment.

           Warning: This function does not consider what happens if a chain in pdb_id_1 matches two chains in pdb_id_2
                    which have differing sequences. In this case, an self-inconsistent set of mutations is returned.
                    One solution would be to extend the mapping to:

                        chain_m_of_pdb_1 -> common_mutations -> List[PDBMutationPair]
                                         -> chain_x_of_pdb_2 -> List[PDBMutationPair]
                                         -> chain_y_of_pdb_2 -> List[PDBMutationPair]
                                         ...

                    where common_mutations contains the set of mutations common to the mapping m->x and m->y etc. whereas the
                    other two mappings contain the set of mutations from m->x or m->y etc. respectively which do not occur in
                    common_mutations. In general, both chain_x_of_pdb_2 and chain_y_of_pdb_2 will be empty excepting the
                    considered case where x and y differ in sequence.

def get_mutations(self, pdb_id_1, pdb_id_2, alignment_tool = 'clustalw', gap_opening_penalty = 0.2):
        '''Returns a mapping chain_of_pdb_1 -> List[PDBMutationPair] representing the mutations needed to transform each chain of pdb_1 into the respective chains of pdb_2.
           This function also sets self.seqres_sequence_maps, a mapping Tuple(pdb_id_1, chain ID in pdb_id_1) -> Tuple(pdb_id_2, chain ID in pdb_id_2) -> a SequenceMap representing the mapping of residues between the chains based on the alignment.

           Warning: This function does not consider what happens if a chain in pdb_id_1 matches two chains in pdb_id_2
                    which have differing sequences. In this case, an self-inconsistent set of mutations is returned.
                    One solution would be to extend the mapping to:

                        chain_m_of_pdb_1 -> common_mutations -> List[PDBMutationPair]
                                         -> chain_x_of_pdb_2 -> List[PDBMutationPair]
                                         -> chain_y_of_pdb_2 -> List[PDBMutationPair]
                                         ...

                    where common_mutations contains the set of mutations common to the mapping m->x and m->y etc. whereas the
                    other two mappings contain the set of mutations from m->x or m->y etc. respectively which do not occur in
                    common_mutations. In general, both chain_x_of_pdb_2 and chain_y_of_pdb_2 will be empty excepting the
                    considered case where x and y differ in sequence.
        '''

        # Set up the objects
        p1 = self.add_pdb(pdb_id_1)
        p2 = self.add_pdb(pdb_id_2)
        self.chain_map[pdb_id_1] = self.chain_map.get(pdb_id_1, {})

        # Determine which chains map to which
        alignment = self.get_alignment(pdb_id_1, pdb_id_2, alignment_tool = alignment_tool, gap_opening_penalty = gap_opening_penalty)
        best_matches = alignment['best_matches']

        # Create the list of mutations
        mutations = {}
        for from_chain, mtches in sorted(best_matches.iteritems()):

            from_pdb_id, from_chain_id = from_chain.split('_')

            # Only consider matches from pdb_id_1 to pdb_id_2
            if from_pdb_id == pdb_id_1:

                self.seqres_sequence_maps[(from_pdb_id, from_chain_id)] = {}

                self.chain_map[from_pdb_id][from_chain_id] = self.chain_map[from_pdb_id].get(from_chain_id, {})

                # Do not consider matches from pdb_id_1 to itself or matches with poor sequence similarity
                restricted_mtchs = {}
                for to_chain, similarity in sorted(mtches.iteritems()):
                    if to_chain.split('_')[0] == pdb_id_2 and similarity >= self.acceptable_sequence_similarity:
                        restricted_mtchs[to_chain] = similarity

                # Take the best matching chains and create a list of mutations needed to transform from_chain to those chains
                # Warning: This does NOT take into account whether the sequences of the best matches differ.
                if restricted_mtchs:
                    top_similarity = max(restricted_mtchs.values())

                    #todo: if the sequences of the best matches differ, raise an Exception. Use 2ZNW and 1DQJ as an example (2ZNW chain A matches with 48% to both 1DQJ chain A and chain B)
                    #top_matches = [to_chain for to_chain, similarity in sorted(restricted_mtchs.iteritems()) if similarity == top_similarity]
                    #pprint.pprint(restricted_mtchs)
                    #print(from_pdb_id, from_chain, 'top_matches', top_matches)
                    #sys.exit(0)


                    for to_chain, similarity in sorted(restricted_mtchs.iteritems()):
                        to_pdb_id, to_chain_id = to_chain.split('_')
                        if similarity == top_similarity:
                            #print(from_pdb_id, from_chain_id)
                            #print(restricted_mtchs)
                            #print(to_pdb_id, to_chain, similarity)
                            self.chain_map[from_pdb_id][from_chain_id][to_pdb_id] = self.chain_map[from_pdb_id][from_chain_id].get(to_pdb_id, set())
                            self.chain_map[from_pdb_id][from_chain_id][to_pdb_id].add(to_chain)
                            mutations[from_chain_id] = mutations.get(from_chain_id, [])
                            chain_mutations = self.get_chain_mutations(from_pdb_id, from_chain_id, to_pdb_id, to_chain_id)
                            mutations[from_chain_id].extend(chain_mutations)

        # mutations can contain duplicates so we remove those
        for chain_id, mlist in mutations.iteritems():
            mutations[chain_id] = sorted(set(mlist))
        return mutations