This function iterates through each filepath in fasta_file_contents and returns a dict mapping (pdbid, chain, file_name) tuples to sequences:
          - options is the OptionParser member;
          - fasta_file_contents is a map from input filenames to the associated FASTA file contents.

def parse_FASTA_files(options, fasta_file_contents):
    ''' This function iterates through each filepath in fasta_file_contents and returns a dict mapping (pdbid, chain, file_name) tuples to sequences:
          - options is the OptionParser member;
          - fasta_file_contents is a map from input filenames to the associated FASTA file contents.
    '''
    records = {}
    reverse_mapping = {}
    original_segment_list = []
    key_location = {}
    sequenceLine = re.compile("^[A-Z]+\n?$")
    sequence_offsets = {}

    for fasta_file_name, tagged_fasta in sorted(fasta_file_contents.iteritems()):

        # Check the tagged pair
        fasta = tagged_fasta[0].strip().split('\n')
        file_type = tagged_fasta[1]
        assert(file_type == 'PDB' or file_type == 'FASTA')
        if not fasta:
            raise Exception("Empty FASTA file.")

        first_line = [line for line in fasta if line.strip()][0]
        if first_line[0] != '>':
            raise Exception("The FASTA file %s is not formatted properly - the first non-blank line is not a description line (does not start with '>')." % fasta_file_name)

        key = None
        line_count = 0
        record_count = 0
        file_keys = []
        unique_keys = {}
        for line in fasta:
            line_count += 1
            line = line.strip()
            if line:
                if line[0] == '>':
                    record_count += 1
                    tokens = [t.strip() for t in line[1:].split('|') if t.strip()]
                    if len(tokens) < 2:
                        raise Exception("The description line ('%s') of record %d of %s is invalid. It must contain both a protein description and a chain identifier, separated by a pipe ('|') symbol." % (line, record_count, fasta_file_name))
                    if len(tokens[0]) < 4:
                        raise Exception("The protein description in the description line ('%s') of record %d of %s is too short. It must be at least four characters long." % (line, record_count, fasta_file_name))
                    if len(tokens[1]) != 1:
                        raise Exception("The chain identifier in the description line ('%s') of record %d of %s is the wrong length. It must be exactky one character long." % (line, record_count, fasta_file_name))

                    # Note: We store the PDB ID as lower-case so that the user does not have to worry about case-sensitivity here (we convert the user's PDB ID argument to lower-case as well)
                    key = (tokens[0][0:4].lower(), tokens[1], fasta_file_name)
                    sub_key = (key[0], key[1]) # this is the part of the key that we expect to be unique (the actual key)
                    key_location[key] = fasta_file_name
                    if sub_key in unique_keys:
                        # todo: we include the fasta_file_name in the key - should we not be checking for uniqueness w.r.t. just tokens[0][0:4].lower() and tokens[1] i.e. omitting the fasta_file_name as part of the check for a more stringent check?
                        raise Exception("Duplicate protein/chain identifier pair. The key %s was generated from both %s and %s. Remember that the first four characters of the protein description are concatenated with the chain letter to generate a 5-character ID which must be unique." % (key, key_location[key], fasta_file_name))
                    records[key] = [line]
                    unique_keys[sub_key] = True
                    file_keys.append(key)
                else:
                    mtchs = sequenceLine.match(line)
                    if not mtchs:
                        raise FastaException("Expected a record header or sequence line at line %d." % line_count)
                    records[key].append(line)

        offset = 0
        if file_type == 'PDB':
            for key in file_keys:
                sequence_length = len(''.join(records[key][1:]))
                sequence_offsets[key[0] + key[1]] = (offset, offset + 1, offset + sequence_length) # storing the sequence start and end residue IDs here is redundant but simplifies code later on
                offset += sequence_length

    # We remove non-protein chains from fragment generation although we did consider them above when determining the offsets
    # as we expect them to be used in predictions
    non_protein_records = []
    set_of_rna_dna_codes = set(('A', 'C', 'G', 'T', 'U', 'X', 'Z'))
    for key, content_lines in records.iteritems():
        mm_sequence = ''.join(content_lines[1:])
        assert(re.match('^[A-Z]+$', mm_sequence)) # Allow X or Z because these may exist (X from the RCSB, Z from our input files)
        if set(mm_sequence).union(set_of_rna_dna_codes) == set_of_rna_dna_codes:
            non_protein_records.append(key)
    for non_protein_record in non_protein_records:
        del records[non_protein_record]

    # If a loops file was passed in, use that to cut up the sequences and concatenate these subsequences to generate a
    # shorter sequence to process. This should save a lot of time when the total length of the subsequences is considerably
    # shorter than the length of the total sequence e.g. in cases where the protein has @1000 residues but we only care about
    # 100 residues in particular loop regions.

    # We need to sample all sequences around a loop i.e. if a sequence segment is 7 residues long at positions 13-19 and we
    # require 9-mers, we must consider the segment from positions 5-27 so that all possible 9-mers are considered.
    residue_offset = max(options.frag_sizes)

    if options.loops_file or options.indices:
        loops_definition = None
        if options.loops_file:
            loops_definition = LoopsFile.from_filepath(options.loops_file, ignore_whitespace = True, ignore_errors = False)

        # If the user supplied more ranges of residues, use those as well
        if options.indices:
            if not loops_definition:
                loops_definition = LoopsFile('')
            for p in options.indices:
                if loops_definition:
                    loops_definition.add(p[0], p[1])

        segment_list = loops_definition.get_distinct_segments(residue_offset, residue_offset)
        original_segment_list = loops_definition.get_distinct_segments(1, 1) # We are looking for 1-mers so the offset is 1 rather than 0

        # Sanity checks
        assert(sorted(segment_list) == segment_list) # sanity check
        for x in range(len(segment_list)):
            segment = segment_list[x]
            if x < len(segment_list) - 1:
                assert(segment[1] < segment_list[x+1][0]) # sanity check

        # Create the generic reverse_mapping from the indices in the sequences defined by the segment_list to the indices in the original sequences.
        # This will be used in FASTA sequences to rewrite the fragments files to make them compatible with the original sequences.
        # Note that this mapping ignores the length of the sequences (which may vary in length) so it may be mapping residues indices
        # which are outside of the length of some of the sequences.
        # Create a sorted list of residues of the chain that we will be including in the sequence for fragment generation
        # then turn that into a 1-indexed mapping from the order of the residue in the sequence to the original residue ID in the PDB
        residues_for_generation = []
        for s in segment_list:
            residues_for_generation += range(s[0], s[1] + 1)
        reverse_mapping['FASTA'] = dict((key, value) for (key, value) in zip(range(1, len(residues_for_generation) + 1), residues_for_generation))

        # Create the reverse_mappings from the indices in the PDB sequences defined by the segment_list to the indices in the original sequences.
        # Membership in sequence_offsets implies a PDB sequence.
        for k, v in sorted(sequence_offsets.iteritems()):
            # For each PDB chain, we consider the set of segments (ignoring extra residues due to nmerage for now so that
            # we do not include chains by accident e.g. if the user specified the first residues of chain C but none in chain B,
            # they probably do not wish to generate fragments for chain B)
            chain_residues = range(v[1], v[2] + 1)
            residues_for_generation = []
            for s in original_segment_list:
                # If the original segment lists lie inside the chain residues then we extend the range w.r.t. the nmers
                if (chain_residues[0] <= s[0] <= chain_residues[-1]) or (chain_residues[0] <= s[1] <= chain_residues[-1]):
                    residues_for_generation += range(s[0] - residue_offset + 1, s[1] + residue_offset - 1 + 1)

            # Create a sorted list of residues of the chain that we will be including in the sequence for fragment generation
            # then turn that into a 1-indexed mapping from the order of the residue in the sequence to the original residue ID in the PDB
            chain_residues_for_generation = sorted(set(chain_residues).intersection(set(residues_for_generation)))
            reverse_mapping[k] = dict((key, value) for (key, value) in zip(range(1, len(chain_residues_for_generation) + 1), chain_residues_for_generation))

        found_at_least_one_sequence = False
        for k, v in sorted(records.iteritems()):
            assert(v[0].startswith('>'))
            subkey = k[0] + k[1]
            sequence = ''.join([s.strip() for s in v[1:]])
            assert(sequenceLine.match(sequence) != None) # sanity check
            cropped_sequence = None
            if sequence_offsets.get(subkey):
                # PDB chain case
                first_residue_id = sequence_offsets[subkey][1]
                cropped_sequence = ''.join([sequence[rmv - first_residue_id] for rmk, rmv in sorted(reverse_mapping[subkey].iteritems())])
                # Sanity check - check that the remapping from the cropped sequence to the original sequence will work in postprocessing
                for x in range(0, len(cropped_sequence)):
                    assert(cropped_sequence[x] == sequence[reverse_mapping[subkey][x + 1] - sequence_offsets[subkey][0] - 1])
                records[k] = [v[0]] + [cropped_sequence[i:i+60] for i in range(0, len(cropped_sequence), 60)] # update the record to only use the truncated sequence
            else:
                # FASTA chain case
                cropped_sequence = ''.join([sequence[rmv - 1] for rmk, rmv in sorted(reverse_mapping['FASTA'].iteritems()) if rmv <= len(sequence)])
                # Sanity check - check that the remapping from the cropped sequence to the original sequence will work in postprocessing
                for x in range(0, len(cropped_sequence)):
                    assert(cropped_sequence[x] == sequence[reverse_mapping['FASTA'][x + 1] - 1])

            found_at_least_one_sequence = found_at_least_one_sequence or (not not cropped_sequence)
            if cropped_sequence:
                records[k] = [v[0]] + [cropped_sequence[i:i+60] for i in range(0, len(cropped_sequence), 60)]
            else:
                del records[k] # delete the chain. todo: test that this works

        if not found_at_least_one_sequence:
            raise Exception('No sequences were created from the loops/indices and the input sequences. This may be an input error so the job is being terminated.')

    if reverse_mapping:
        return records, dict(reverse_mapping = reverse_mapping, segment_list = original_segment_list, sequence_offsets = sequence_offsets)
    else:
        return records, None