docs: 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