def METHOD_NAME():
    parser = argparse.ArgumentParser( description='Report coverage gaps/dips from a samtools mpileup file')
    ## output file to be written
    parser.add_argument('-i', '--input_file', type=str, required=True, help='Path to an input mpileup file to be read' )
    parser.add_argument('-o', '--output_file', type=str, required=False, help='Path to an output file to be created' )
    parser.add_argument('-f', '--fasta_file', type=str, required=True, help='Reference fasta file, against which reads were aligned.  Needed for low 3-prime end coverage' )
    parser.add_argument('-mcd', '--min_coverage_depth', type=int, required=True, help='Min coverage depth, below which is reported' )
    parser.add_argument('-mcs', '--min_coverage_span', type=int, required=False, help='Coverage window size, the avg of which is calculated for depth cutoff' )
    parser.add_argument('-eb', '--end_buffer', type=int, required=False, default=0, help='If specified, gaps this far from either end of the molecule will not be reported.' )
    args = parser.parse_args()
    # Check, this isn't ready yet:
    if args.min_coverage_span is not None:
        raise Exception("ERROR: Sorry, --min_coverage_span not yet implemented.")
    if args.output_file is None:
        out_fh = sys.stdout
    else:
        out_fh = open(args.output_file, 'wt')
    lengths = utils.fasta_sizes_from_file(args.fasta_file)
    stats = {'total_molecules': 0, 'total_bases': 0, 'depth_sum': 0}
    # In mpileup gaps are reported either with positions of coverage 0 OR omitted rows.
    depths = list()
    current_seq_id = None
    for line in open(args.input_file):
        contig, this_coord, base, depth = line.split("\t")[0:4]
        this_coord = int(this_coord)
        # new molecule
        if contig != current_seq_id:
            stats['total_molecules'] += 1
            # purge the last one
            if current_seq_id != None:
                print_spans(current_seq_id, depths, args.min_coverage_depth, out_fh, stats, lengths, args.end_buffer)
            depths = [0] * lengths[contig]
            current_seq_id = contig
        depths[this_coord - 1] = depth
    print_spans(current_seq_id, depths, args.min_coverage_depth, out_fh, stats, lengths, args.end_buffer)
    print("INFO: Total molecules: {0}".format(stats['total_molecules']), file=sys.stderr)
    print("INFO: Total bases    : {0}".format(stats['total_bases']), file=sys.stderr)
    print("INFO: Avg cov depth  : {0}x".format(int(stats['depth_sum'] / stats['total_bases'])), file=sys.stderr)