text: """Contains functions to parse and preprocess information from the input file"""
  import sys
  import os
  import h5py
  import logging
  import multiprocessing as mp
  import numpy as np
  import pickle
  import signal as sig
  
  from immunopepper.io_ import decodeUTF8
  from immunopepper.namedtuples import CountInfo
  from immunopepper.namedtuples import GeneInfo
  from immunopepper.namedtuples import GeneTable
  from immunopepper.namedtuples import ReadingFrameTuple
  from immunopepper.utils import encode_chromosome
  from immunopepper.utils import find_overlapping_cds_simple
  from immunopepper.utils import get_successor_list
  from immunopepper.utils import leq_strand
  from immunopepper.mutations import get_mutation_mode_from_parser
  from immunopepper.mutations import parse_mutation_file
  
  def genes_preprocess_batch(genes, gene_idxs, gene_cds_begin_dict, verbose=False):
      gene_info = []
      for gene in genes:
          gene.from_sparse()
          assert (gene.strand in ["+", "-"])
          assert (len(gene.transcripts) == len(gene.exons))
  
          # Ignore genes that have no CDS annotated...
          if gene.name not in gene_cds_begin_dict:
              gene_info.append(None)
              continue
  
          vertex_succ_list = get_successor_list(gene.splicegraph.edges, gene.splicegraph.vertices, gene.strand)
          if gene.strand == "+":
              vertex_order = np.argsort(gene.splicegraph.vertices[0, :])
          else:  # gene.strand=="-"
              vertex_order = np.argsort(gene.splicegraph.vertices[1, :])[::-1]
  
          # get the reading_frames
          reading_frames = {}
          vertex_len_dict = {}
          for idx in vertex_order:
              reading_frames[idx] = set()
              v_start = gene.splicegraph.vertices[0, idx]
              v_stop = gene.splicegraph.vertices[1, idx]
              cds_begins = find_overlapping_cds_simple(v_start, v_stop, gene_cds_begin_dict[gene.name], gene.strand)
              vertex_len_dict[idx] = v_stop - v_start
  
              # Initialize reading regions from the CDS transcript annotations
              for cds_begin in cds_begins:
                  line_elems = cds_begin[2]
                  cds_strand = line_elems[6]
                  assert (cds_strand == gene.strand)
                  cds_phase = int(line_elems[7])
                  cds_left = int(line_elems[3]) - 1
                  cds_right = int(line_elems[4])
  
                  # TODO: need to remove the redundance of (cds_start, cds_stop, item)
                  if gene.strand == "-":
                      cds_right_modi = max(cds_right - cds_phase, v_start)
                      cds_left_modi = v_start
                      n_trailing_bases = cds_right_modi - cds_left_modi
                  else:
                      cds_left_modi = min(cds_left + cds_phase, v_stop)
                      cds_right_modi = v_stop
                      n_trailing_bases = cds_right_modi - cds_left_modi
  
                  read_phase = n_trailing_bases % 3
                  reading_frames[idx].add(ReadingFrameTuple(cds_left_modi, cds_right_modi, read_phase))
          gene.to_sparse()
          gene_info.append(GeneInfo(vertex_succ_list, vertex_order, reading_frames, vertex_len_dict,
                                    gene.splicegraph.vertices.shape[1]))
  
      return gene_info, gene_idxs
  
  
  def genes_preprocess_all(genes, gene_cds_begin_dict, parallel=1):
      """ Preprocess the gene and generate new attributes under gene object
          Modify the gene object directly
  
      Parameters
      ----------
      genes: List[Object]. List of gene objects. The object is generated by SplAdder
      gene_cds_begin_dict: Dict. str -> List(int) From gene name to list of cds start positions
      """
  
      if parallel > 1:
          global genes_info
          global cnt
          genes_info = np.zeros((genes.shape[0],), dtype=object)
          cnt = 0
  
          def update_gene_info(result):
              global genes_info
              global cnt
              for i, tmp in enumerate(result[0]):
                  if cnt > 0 and cnt % 100 == 0:
                      sys.stdout.write('.')
                      if cnt % 1000 == 0:
                          sys.stdout.write('%i/%i\n' % (cnt, genes.shape[0]))
                      sys.stdout.flush()
                  cnt += 1
                  genes_info[result[1][i]] = tmp
              del result
  
          pool = mp.Pool(processes=parallel, initializer=lambda: sig.signal(sig.SIGINT, sig.SIG_IGN))
          for i in range(0, genes.shape[0], 100):
              gene_idx = np.arange(i, min(i + 100, genes.shape[0]))
              _ = pool.apply_async(genes_preprocess_batch, args=(genes[gene_idx], gene_idx, gene_cds_begin_dict,),
                                   callback=update_gene_info)
          pool.close()
          pool.join()
      else:
          genes_info = genes_preprocess_batch(genes, np.arange(genes.shape[0]), gene_cds_begin_dict, verbose=True)[0]
  
      return genes_info
  
  
  def preprocess_ann(ann_path):
      """ Extract information from annotation file (.gtf, .gff and .gff3)
  
      Parameters
      ----------
      ann_path: str. Annotation file path
  
      Returns
      -------
      gene_table: NamedTuple.store the gene-transcript-cds mapping tables derived
          from .gtf file. has attribute ['gene_to_cds_begin', 'ts_to_cds', 'gene_to_cds']
      chromosome_set: set. Store the chromosome naming.
      """
      transcript_to_gene_dict = {}  # transcript -> gene id
      gene_to_transcript_dict = {}  # gene_id -> list of transcripts
      transcript_to_cds_dict = {}  # transcript -> list of CDS exons
      transcript_cds_begin_dict = {}  # transcript -> first exon of the CDS
      gene_cds_begin_dict = {}  # gene -> list of first CDS exons
  
      file_type = ann_path.split('.')[-1]
      chromesome_set = set()
      # collect information from annotation file
      for line in open(ann_path, 'r'):
          if line[0] == '#':
              continue
          item = line.strip().split('\t')
          chromesome_set.add(item[0])
          feature_type = item[2]
          attribute_item = item[-1]
          attribute_dict = attribute_item_to_dict(attribute_item, file_type, feature_type)
          # store relationship between gene ID and its transcript IDs
          if feature_type in ['transcript', 'mRNA']:
              gene_id = attribute_dict['gene_id']
              transcript_id = attribute_dict['transcript_id']
              if attribute_dict['gene_type'] != 'protein_coding' or attribute_dict['transcript_type'] != 'protein_coding':
                  continue
              assert (transcript_id not in transcript_to_gene_dict)
              transcript_to_gene_dict[transcript_id] = gene_id
              if gene_id in gene_to_transcript_dict and transcript_id not in gene_to_transcript_dict[gene_id]:
                  gene_to_transcript_dict[gene_id].append(transcript_id)
              else:
                  gene_to_transcript_dict[gene_id] = [transcript_id]
  
          # Todo python is 0-based while gene annotation file(.gtf, .vcf, .maf) is one based
          elif feature_type == "CDS":
              parent_ts = attribute_dict['transcript_id']
              strand_mode = item[6]
              cds_left = int(item[3]) - 1
              cds_right = int(item[4])
              frameshift = int(item[7])
              if parent_ts in transcript_to_cds_dict:
                  transcript_to_cds_dict[parent_ts].append((cds_left, cds_right, frameshift))
              else:
                  transcript_to_cds_dict[parent_ts] = [(cds_left, cds_right, frameshift)]
              if strand_mode == "+":
                  cds_start, cds_stop = cds_left, cds_right
              else:
                  cds_start, cds_stop = cds_right, cds_left
  
              # we only consider the start of the whole CoDing Segment
              if parent_ts not in transcript_cds_begin_dict or \
                      leq_strand(cds_start, transcript_cds_begin_dict[parent_ts][0], strand_mode):
                  transcript_cds_begin_dict[parent_ts] = (cds_start, cds_stop, item)
  
      # collect first CDS exons for all transcripts of a gene
      for ts_key in transcript_to_gene_dict:
          target_gene = transcript_to_gene_dict[ts_key]
          if target_gene not in gene_cds_begin_dict:
              gene_cds_begin_dict[target_gene] = []
          if ts_key in transcript_cds_begin_dict:
              gene_cds_begin_dict[target_gene].append(transcript_cds_begin_dict[ts_key])
  
      # sort list of CDS exons per transcript
      for ts_key in transcript_to_cds_dict:
          transcript_to_cds_dict[ts_key] = sorted(transcript_to_cds_dict[ts_key], key=lambda coordpair: coordpair[0])
  
      genetable = GeneTable(gene_cds_begin_dict, transcript_to_cds_dict, gene_to_transcript_dict)
      return genetable, chromesome_set
  
  
  def attribute_item_to_dict(a_item, file_type, feature_type):
      """  From attribute item in annotation file to get corresponding dictionary
  
      Parameters
      ----------
      a_item: str. attribute item
      file_type: str. Choose from {'gtf', 'gff', 'gff3'}
      feature_type: str. Extract other fields. We only
          consider 'CDS', 'mRNA' and 'transcript'
  
      Returns
      -------
      gtf_dict: dict. store all the necessary data
  
      """
      gtf_dict = {}
      if file_type.lower() == 'gtf':
          attribute_list = a_item.split('; ')
          for attribute_pair in attribute_list:
              pair = attribute_pair.split(' ')
              gtf_dict[pair[0]] = pair[1][1:-1]
      elif file_type.lower() == 'gff3':
          attribute_list = a_item.split(';')
          for attribute_pair in attribute_list:
              pair = attribute_pair.split('=')
              gtf_dict[pair[0]] = pair[1]
      elif file_type.lower() == 'gff':
          gff_dict = {}
          attribute_list = a_item.split(';')
          for attribute_pair in attribute_list:
              pair = attribute_pair.split('=')
              gff_dict[pair[0]] = pair[1]  # delete "", currently now work on level 2
          if feature_type == 'CDS':
              gtf_dict['transcript_id'] = gff_dict['Parent']
          elif feature_type in {'mRNA', 'transcript'}:  # mRNA or transcript
              gtf_dict['gene_id'] = gff_dict['geneID']
              gtf_dict['transcript_id'] = gff_dict['ID']
              gtf_dict['gene_type'] = gff_dict['gene_type']
              gtf_dict['transcript_type'] = gff_dict['transcript_type']
  
      return gtf_dict
  
  
  def search_edge_metadata_segmentgraph(gene, coord, countinfo, Idx, edge_idxs=None, edge_counts=None):
      """Given the ordered edge coordinates of the edge, return expression information of the edge
  
      Parameters
      ----------
      gene: Object. Generated by SplAdder
      coord: np.array of length 4. Sorted coordinates of 4 positions in ascending order
      countinfo: NamedTuple, contains SplAdder count info
      Idx: Namedtuple, has attribute idx.gene and idx.sample
      edge_idxs: np.array, containing the edge index values for the current gene
      egde_counts: np.array, containing the edge count values for the current gene
  
      Returns
      -------
      count: tuple of float. Expression level for the given edges.
      """
  
      def get_segmentgraph_edge_expr(sorted_pos, h5f, edge_idxs, edge_counts=None):
          a = np.searchsorted(segmentgraph.segments[1, :], sorted_pos[1])
          b = np.searchsorted(segmentgraph.segments[0, :], sorted_pos[2])
          if a < b:
              idx = np.ravel_multi_index([a, b], segmentgraph.seg_edges.shape)
          else:
              idx = np.ravel_multi_index([b, a], segmentgraph.seg_edges.shape)
          cidx = np.searchsorted(edge_idxs, idx)
          if edge_counts is None:
              return h5f['edges'][edge_gene_idxs[cidx], Idx.sample]
          else:
              return edge_counts[cidx]
  
      segmentgraph = gene.segmentgraph
      sorted_pos = np.sort(np.array([coord.start_v1, coord.stop_v1, coord.start_v2, coord.stop_v2]))
  
      h5f = None
      if edge_counts is None:
          h5f = h5py.File(countinfo.h5fname, 'r')
  
      if edge_idxs is None or edge_counts is None:
          edge_gene_idxs = np.where(countinfo.gene_ids_edges == countinfo.gene_idx_dict[gene.name])[0]
          edge_idxs = h5f['edge_idx'][edge_gene_idxs].astype('int')
          count = get_segmentgraph_edge_expr(sorted_pos, h5f, edge_idxs)
      else:
          count = get_segmentgraph_edge_expr(sorted_pos, h5f, edge_idxs, edge_counts)
  
      if coord.start_v3 is None:
          return (count,)
      else:
          sorted_pos = np.sort(np.array([coord.start_v2, coord.stop_v2, coord.start_v3, coord.stop_v3]))
          if edge_idxs is None or edge_counts is None:
              count2 = get_segmentgraph_edge_expr(sorted_pos, h5f, edge_idxs)
          else:
              count2 = get_segmentgraph_edge_expr(sorted_pos, h5f, edge_idxs, edge_counts)
          return (count, count2)
  
  
  def parse_gene_metadata_info(h5fname, sample_list):
      """ Parse the count file
  
      Parameters
      ----------
      h5fname: str. .count.h5f file
      sample_list: List(str). List of samples.
  
      Returns
      -------
      countinfo: Namedtuple. Store all the counts information. Has attributes:
          'sample_idx_dict' --> dict from sample name to index
          'gene_idx_dict' --> dict from gene name to index
          'gene_ids_segs' --> array containing segment-geneID relation
          'gene_ids_edges' --> array containing edge-geneID relation
          'h5fname' --> HDF5 file name
      """
  
      # the SplAdder count hdf5 file has the following structure
      #   h5f["strains"] --> sample names
      #   h5f["segments"] --> segment expression (rows: segments, columns: samples)
      #   h5f["edges"] --> edge expression (rows: edges, columns: samples)
      #   h5f["edge_idx"] --> multi-row index encoding the edge in the splice graph (rows:
      h5f = h5py.File(h5fname, 'r')
      assert (h5f["strains"].shape[0] == h5f["segments"].shape[1])
      assert (h5f["gene_ids_segs"].size == h5f["segments"].shape[0])
      assert (h5f["gene_ids_edges"].size == h5f["edges"].shape[0])
  
      ### create a sample name dictionary mapping sample names to indices
      sample_names = h5f['strains'][:] if len(h5f['strains'].shape) == 1 else h5f['strains'][:, 0]
      sample_idx_dict = dict([(n.decode('utf8'), i) for i, n in enumerate(sample_names)])
  
      ### create a gene name dictionary mapping gene names to indices
      gene_names = h5f['gene_names'][:] if len(h5f['gene_names'].shape) == 1 else h5f['gene_names'][:, 0]
      gene_idx_dict = dict([(n.decode('utf8'), i) for i, n in enumerate(gene_names)])
  
      gene_ids_segs = h5f['gene_ids_segs'][:] if len(h5f['gene_ids_segs'].shape) == 1 else h5f['gene_ids_segs'][:, 0]
      gene_ids_edges = h5f['gene_ids_edges'][:] if len(h5f['gene_ids_edges'].shape) == 1 else h5f['gene_ids_edges'][:, 0]
  
      ### segs
      gene_ids_segs_u, gene_ids_segs_idx = np.unique(gene_ids_segs, return_index=True)
      gene_ids_segs_idx_last = np.r_[gene_ids_segs_idx[1:], [gene_ids_segs.shape[0]]]
      gene_id_to_segrange = dict()
      for i, g in enumerate(gene_ids_segs_u):
          gene_id_to_segrange[g] = (gene_ids_segs_idx[i], gene_ids_segs_idx_last[i])
      ### edges
      gene_ids_edges_u, gene_ids_edges_idx = np.unique(gene_ids_edges, return_index=True)
      gene_ids_edges_idx_last = np.r_[gene_ids_edges_idx[1:], [gene_ids_edges.shape[0]]]
      gene_id_to_edgerange = dict()
      for i, g in enumerate(gene_ids_edges_u):
          gene_id_to_edgerange[g] = (gene_ids_edges_idx[i], gene_ids_edges_idx_last[i])
  
      countinfo = CountInfo(sample_idx_dict,
                            gene_idx_dict,
                            gene_id_to_segrange,
                            gene_id_to_edgerange,
                            h5fname)
      h5f.close()
      return countinfo
  
  
  def parse_mutation_from_vcf(vcf_path, output_dir='', heter_code=0, mut_pickle=False, h5_sample_list=None):
      """Extract germline mutation information from the given vcf file and vcf.h5 file
  
      Parameters
      ----------
      vcf_path: str, vcf file path
      output_dir: str, path to vcf pickle output directory
      heter_code: int (0 or 2). specify which number represents heter alle.
          0: 0-> homozygous alternative(1|1), 1-> heterozygous(0|1,1|0) 2->homozygous reference(0|0)
          2: 0-> homozygous reference(0|0), 1-> heterozygous(0|1,1|0) 2->homozygous alternative(1|1)
      mut_pickle: bool, flag indicating whether to pickle mutation info to disk
      h5_sample_list: list, list of samples
  
      Returns
      -------
      mut_dict: with key (sample, chromo) and values (var_dict)
      """
      vcf_pkl_file = os.path.join(output_dir, 'vcf.pickle')
      if mut_pickle:
          if os.path.exists(vcf_pkl_file):
              f = open(vcf_pkl_file, 'rb')
              mutation_dic = pickle.load(f)
              logging.info("Use pickled vcf mutation dict in {}".format(vcf_pkl_file))
              return mutation_dic
  
      file_type = vcf_path.split('.')[-1]
      if file_type == 'h5':  # hdf5 filr
          mutation_dic = parse_mutation_from_vcf_h5(vcf_path, h5_sample_list, heter_code)
          logging.info("Get germline mutation dict from h5 file in {}. No pickle file created".format(vcf_path))
          return mutation_dic
      else:  # vcf text file
          f = open(vcf_path, 'r')
          lines = f.readlines()
          mutation_dic = {}
          for line in lines:
              if line.strip()[:2] == '##':  # annotation line
                  continue
              if line.strip()[0] == '#':  # head line
                  fields = line.strip().split('\t')
                  sample_list = fields[9:]
                  continue
              items = line.strip().split('\t')
              var_dict = {}
              chr = items[0]
              pos = int(items[1]) - 1
              var_dict['ref_base'] = items[3]
              var_dict['mut_base'] = items[4]
              var_dict['qual'] = items[5]
              var_dict['filter'] = items[6]
              if len(var_dict['ref_base']) == len(var_dict['mut_base']):  # only consider snp for now
                  for i, sample_id in enumerate(sample_list):
                      if items[9 + i].split(':')[0] in {'1|1', '1|0', '0|1'}:
                          if (sample_id, chr) in list(mutation_dic.keys()):
                              mutation_dic[(sample_id, chr)][int(pos)] = var_dict
                          else:
                              mutation_dic[(sample_id, chr)] = {}
                              mutation_dic[(sample_id, chr)][int(pos)] = var_dict
      if mut_pickle:
          f_pkl = open(vcf_pkl_file, 'wb')
          pickle.dump(mutation_dic, f_pkl)
          logging.info("create vcf pickled mutation dict for next time's use in {}".format(vcf_pkl_file))
  
      return mutation_dic
  
  
  def parse_mutation_from_vcf_h5(h5_vcf_path, sample_list, heter_code=0):
      """
      Extract germline mutation information from given vcf h5py file.
  
      Parameters
      ----------
      h5_vcf_path: str, vcf file path
      sample_list: list of str, list for sample name
      heter_code: int (0 or 2). specify which number represents heter alle.
          0: 0-> homozygous alternative(1|1), 1-> heterozygous(0|1,1|0) 2->homozygous reference(0|0)
          2: 0-> homozygous reference(0|0), 1-> heterozygous(0|1,1|0) 2->homozygous alternative(1|1)
  
      Returns
      -------
      mut_dict: with key (sample, chromo) and values (var_dict)
  
      """
      a = h5py.File(h5_vcf_path, 'r')
      mut_dict = {}
      for sample in sample_list:
          col_id = [i for (i, item) in enumerate(a['gtid']) if decodeUTF8(item).startswith(sample)][0]
          row_id = np.where(np.logical_or(a['gt'][:, col_id] == heter_code, a['gt'][:, col_id] == 1))[0]
          for irow in row_id:
              chromo = encode_chromosome(a['pos'][irow, 0])
              pos = a['pos'][irow, 1] - 1
              mut_base = decodeUTF8(a['allele_alt'][irow])
              ref_base = decodeUTF8(a['allele_ref'][irow])
              var_dict = {"mut_base": mut_base, "ref_base": ref_base}
              if (sample, chromo) in mut_dict:
                  mut_dict[(sample, chromo)][pos] = var_dict
              else:
                  mut_dict[(sample, chromo)] = {}
                  mut_dict[(sample, chromo)][pos] = var_dict
      return mut_dict
  
  
  def parse_mutation_from_maf(maf_path, output_dir='', mut_pickle=False):
      """
      Extract somatic mutation information from given maf file.
  
      Parameters
      ----------
      maf_path: str, maf file path
      output_dir: str, save a pickle for maf_dict to save preprocess time
      mut_pickle: bool, flag indicating whether to pickle mutation info to disk
  
      Returns
      -------
      mut_dict: with key (sample, chromo) and values (var_dict)
  
      """
      maf_pkl_file = os.path.join(output_dir, 'maf.pickle')
      if mut_pickle:
          if os.path.exists(maf_pkl_file):
              f = open(maf_pkl_file, 'rb')
              mutation_dic = pickle.load(f)
              logging.info("Use pickled maf mutation dict in {}".format(maf_pkl_file))
              return mutation_dic
  
      f = open(maf_path)
      lines = f.readlines()
      mutation_dic = {}
      for i, line in enumerate(lines[1:]):
          print(i)
          items = line.strip().split('\t')
          if items[9] == 'SNP':  # only consider snp
              sample_id = items[15]
              chr = items[4]
              pos = int(items[5]) - 1
              var_dict = {}
              var_dict['ref_base'] = items[10]
              var_dict['mut_base'] = items[12]
              var_dict['strand'] = items[7]
              var_dict['variant_Classification'] = items[8]
              var_dict['variant_Type'] = items[9]
              if (sample_id, chr) in list(mutation_dic.keys()):
                  mutation_dic[((sample_id, chr))][int(pos)] = var_dict
              else:
                  mutation_dic[((sample_id, chr))] = {}
                  mutation_dic[((sample_id, chr))][int(pos)] = var_dict
      if mut_pickle:
          f_pkl = open(maf_pkl_file, 'wb')
          pickle.dump(mutation_dic, f_pkl)
          logging.info("create maf pickled mutation dict for next time's use in {}".format(maf_pkl_file))
      return mutation_dic
  
  
  # todo: support tsv file in the future
  def parse_junction_meta_info(h5f_path):
      """ Extract introns of interest from given h5py file
  
      Parameters
      ----------
      h5f_path: str, h5py file path
  
      Returns
      -------
      junction_dict: dict, key (chromosome id), value (set of coordinate pairs)
  
      """
      if h5f_path is None:
          return None
      else:
          h5f = h5py.File(h5f_path, 'r')
          chrms = h5f['chrms'][:]
          pos = h5f['pos'][:].astype('str')
          strand = h5f['strand'][:]
          junction_dict = {}
  
          for i, ichr in enumerate(chrms):
              try:
                  junction_dict[decodeUTF8(ichr)].add(':'.join([pos[i, 0], pos[i, 1], decodeUTF8(strand[i])]))
              except KeyError:
                  junction_dict[decodeUTF8(ichr)] = set([':'.join([pos[i, 0], pos[i, 1], decodeUTF8(strand[i])])])
      return junction_dict
  
  
  
  def main():
      base_path = '/Users/laurieprelot/Documents/Projects/immunopepper/data_test'
      mutation_mode = 'somatic'
      germline_file_path = os.path.join(base_path, 'mergedfiles_clean_stringentfilter.matchIds.h5')
      somatic_file_path = os.path.join(base_path, 'pancan.merged.v0.2.6.PUBLIC.matchIds.maf')
      output_dir = base_path
      heter_code = 0
      mut_pickle = False
      h5_sample_list = ['TCGA-13-1489-01A-01']
      is_error = True
      if mutation_mode == 'somatic_and_germline':
          if somatic_file_path != '' and germline_file_path != '':
              somatic_mutation_dict = parse_mutation_file(somatic_file_path,output_dir,heter_code,mut_pickle,h5_sample_list)
              germline_mutation_dict = parse_mutation_file(germline_file_path,output_dir,heter_code,mut_pickle,h5_sample_list)
              is_error = False
      elif mutation_mode == 'germline':
          if germline_file_path != '':
              somatic_mutation_dict = {}  # empty dic
              germline_mutation_dict = parse_mutation_file(germline_file_path,output_dir,heter_code,mut_pickle,h5_sample_list)
              is_error = False
      elif mutation_mode == 'somatic':
          if somatic_file_path != '':
              somatic_mutation_dict = parse_mutation_file(somatic_file_path,output_dir,heter_code,mut_pickle,h5_sample_list)
              germline_mutation_dict = {}
              is_error = False
      elif mutation_mode == 'ref':
          somatic_mutation_dict = {}
          germline_mutation_dict = {}
          is_error = False
  
  if __name__ == "__main__":
      main()