docs: def cluster(args):
      """
      Creating clusters
      """
  
      args = _check_args(args)
      read_stats_file = op.join(args.dir_out, "read_stats.tsv")
      if file_exists(read_stats_file):
          os.remove(read_stats_file)
  
      bam_file, seq_obj = _clean_alignment(args)
  
      logger.info("Parsing matrix file")
      seqL, y, l = parse_ma_file(seq_obj, args.ffile)
      # y, l = _total_counts(seqL.keys(), seqL)
      logger.info("counts after: %s" % sum(y.values()))
      logger.info("# sequences after: %s" % l)
      dt = pd.DataFrame({'sample': y.keys(), 'counts': y.values()})
      dt['step'] = 'aligned'
      dt.to_csv(read_stats_file, sep="\t", index=False, header=False, mode='a')
  
      if len(seqL.keys()) < 10:
          logger.error("It seems you have low coverage. Please check your fastq files have enough sequences.")
          raise ValueError("So few sequences.")
  
      logger.info("Cleaning bam file")
      y, l = _total_counts(seqL.keys(), seqL)
      logger.info("counts after: %s" % sum(y.values()))
      logger.info("# sequences after: %s" % l)
      dt = pd.DataFrame({'sample': y.keys(), 'counts': y.values()})
      dt['step'] = 'cleaned'
      dt.to_csv(read_stats_file, sep="\t", index=False, header=False, mode='a')
  
      clusL = _create_clusters(seqL, bam_file, args)
      y, l = _total_counts(clusL.seq.keys(), clusL.seq, aligned=True)
      logger.info("counts after: %s" % sum(y.values()))
      logger.info("# sequences after: %s" % l)
      dt = pd.DataFrame({'sample': y.keys(), 'counts': y.values()})
      dt['step'] = 'clusters'
      dt.to_csv(read_stats_file, sep="\t", index=False, header=False, mode='a')
  
      logger.info("Solving multi-mapping events in the network of clusters")
      clusLred = _cleaning(clusL, args.dir_out)
      y, l = _total_counts(clusLred.clus, seqL)
      logger.info("counts after: %s" % sum(y.values()))
      logger.info("# sequences after: %s" % l)
      dt = pd.DataFrame({'sample': y.keys(), 'counts': y.values()})
      dt['step'] = 'meta-cluster'
      dt.to_csv(read_stats_file, sep="\t", index=False, header=False, mode='a')
      logger.info("Clusters up to %s" % (len(clusLred.clus.keys())))
  
      if args.show:
          logger.info("Creating sequences alignment to precursor")
          clusLred = show_seq(clusLred, args.index)
  
      clusLred = peak_calling(clusLred)
  
      clusLred = _annotate(args, clusLred)
  
      logger.info("Creating json and count matrix")
      json_file = _create_json(clusLred, args)
  
      logger.info("Output file in: %s" % args.dir_out)
  
      if args.db:
          name = args.db + ".db"
          logger.info("Create database: database/" + name)
          data = load_data(json_file)
          out_dir = op.join(args.dir_out, "database")
          make_database(data, name, out_dir)
      logger.info("Finished")