/*
  bed2bed_display.c

  Convert BED file to BED format suitable for displaying ChIP-seq peaks.

  # Arguments:
  #
  #   BED File
  #   Minimal BED score - 5th field (integer)  [-a min]
  #   Maximal BED score - 5th field (integer)  [-b max]
  #
  # Options:
  #
  #   Set Feature name    [-f fname]
  #   BED output formats: [-o 1|2]  def=1
  #
  #     1- BED track with score and strand fields -use score to visualize the track
  #        the score value will determine the level of gray in which this feature is displayed
  #        (higher numbers = darker gray).
  #
  #        Track Parameters:
  #
  #        ex: track name=STAT1Peaks description="Robertson 2007, STAT1 peaks" visibility=2 useScore=1
  #
  #     2- BedGraph format
  #        The bedGraph format allows display of continuous-valued data in track format.
  #        chromA  chromStartA  chromEndA  dataValueA
  #        chromB  chromStartB  chromEndB  dataValueB
  #
  #        Parameters for bedGraph track definition lines:
  #
  #        track type=bedGraph name=track_label description=center_label visibility=1|2|3 color=r,g,b graphType=bar|points

  Giovanna Ambrosini, EPFL/ISREC, Giovanna.Ambrosini@epfl.ch

  Copyright (c) 2014 EPFL and SIB Swiss Institute of Bioinformatics.

  This program is free software: you can redistribute it and/or modify
  it under the terms of the GNU General Public License as published by
  the Free Software Foundation, either version 3 of the License, or
  (at your option) any later version.

  This program is distributed in the hope that it will be useful,
  but WITHOUT ANY WARRANTY; without even the implied warranty of
  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  GNU General Public License for more details.

  You should have received a copy of the GNU General Public License
  along with this program.  If not, see <http://www.gnu.org/licenses/>.

*/
/*
#define DEBUG
*/
#define _GNU_SOURCE
#include <stdlib.h>
#include <stdio.h>
#include <string.h>
#include <errno.h>
#include <unistd.h>
#include <ctype.h>
#include <getopt.h>
#include <assert.h>
#include <sys/stat.h>
#include "hashtable.h"
#ifdef DEBUG
#include <mcheck.h>
#endif

#include "version.h"

#define BUF_SIZE 8192
#define LINE_SIZE 1024
#define FT_MAX  64
#define CHR_NB 64
#define AC_MAX 64
#define POS_MAX 16
#define SCORE_MAX 12
#define NAME_MAX 128
#define FIELD_MAX 64

typedef struct _options_t {
    int help;
    int debug;
    char *feature;
    int oformat;
    char *autoscale;
    char *always0;
    char *wfunction;
    char *smoothing;
    char *visibility;
    char *trackName;
    char *trackDesc;
    char *trackColor;
} options_t;

static options_t options;

int minScore = -1;
int maxScore = -1;

int rescale_score(int s) {      /* Rescale score range to [167,1000] */
    int a = 167;
    int b = 1000;
    int val;

    if (maxScore - minScore > 0) {
        val = ((b - a) * (s - minScore) / (maxScore - minScore)) + a;
    } else {
        val = s;
    }
    return (val);
}

int process_bed(FILE *input, char *iFile) {
    unsigned long start, end;
    int score;
    int new_score;
    int count = 1;
    char *s, *res, *buf;
    size_t bLen = LINE_SIZE;

    if (options.debug)
        fprintf(stderr, " Processing BED file %s\n", iFile);
    if ((s = malloc(bLen * sizeof(char))) == NULL) {
        perror("process_bed: malloc");
        exit(1);
    }
    /* Print Track Header */
    if (options.oformat == 1) {
        printf
            ("track name=\"%s\" description=\"%s\" visibility=%s useScore=1\n",
             options.trackName, options.trackDesc, options.visibility);
    } else if (options.oformat == 2) {
        printf
            ("track type=bedGraph name=\"%s\" description=\"%s\" visibility=%s color=%s autoScale=%s alwaysZero=%s maxHeightPixels=100:50:20 graphType=bar priority=30 windowingFunction=%s smoothingWindow=%s\n",
             options.trackName, options.trackDesc, options.visibility,
             options.trackColor, options.autoscale, options.always0,
             options.wfunction, options.smoothing);
    }
#ifdef DEBUG
    int c = 1;
#endif
    while ((res = fgets(s, (int) bLen, input)) != NULL) {
        size_t cLen = strlen(s);
        char bed_fld[12][FIELD_MAX];
        char strand = '\0';
        char field[32];
        int i = 0;

        while (cLen + 1 == bLen && s[cLen - 1] != '\n') {
            bLen *= 2;
            if ((s = realloc(s, bLen)) == NULL) {
                perror("process_file: realloc");
                exit(1);
            }
            res = fgets(s + cLen, (int) (bLen - cLen), input);
            cLen = strlen(s);
        }
        if (s[cLen - 1] == '\n')
            s[cLen - 1] = 0;
        else if (s[cLen - 1] == '\n' && s[cLen - 2] == '\r')
            s[cLen - 2] = 0;

        buf = s;
        /* printf ("BED line : %s\n", s); */
        /* Check line */
        /* Get first word/field */
        if (sscanf(s, "%31[^ ]", field) == 1) {
            if ((strcmp(field, "track") == 0)
                || (strcmp(field, "browser") == 0)
                || (strcmp(field, "itemRgb") == 0))
                continue;
        }
        /* Get BED fields */
        /* Chrom NB */
        while (*buf != 0 && !isspace(*buf)) {
            if (i >= CHR_NB) {
                bed_fld[0][i-1] = 0;
                fprintf(stderr,
                        "Chrom NB too long \"%s\" Skipping line...\n",
                        bed_fld[0]);
                break;
            }
            bed_fld[0][i++] = *buf++;
        }
        if (i < AC_MAX)
            bed_fld[0][i] = 0;
        else
            continue;
        while (isspace(*buf))
            buf++;
        /* Start Position */
        i = 0;
        while (isdigit(*buf)) {
            if (i >= POS_MAX) {
                fprintf(stderr, "Start position too large \"%s\" \n", buf);
                exit(1);
            }
            bed_fld[1][i++] = *buf++;
        }
        bed_fld[1][i] = 0;
        start = (unsigned long) atoi(bed_fld[1]);
        while (isspace(*buf))
            buf++;
        /* End Position */
        i = 0;
        while (isdigit(*buf)) {
            if (i >= POS_MAX) {
                fprintf(stderr, "End position too large \"%s\" \n", buf);
                exit(1);
            }
            bed_fld[2][i++] = *buf++;
        }
        bed_fld[2][i] = 0;
        end = (unsigned long) atoi(bed_fld[2]);
        while (isspace(*buf))
            buf++;
        /* Name/Feature */
        i = 0;
        while (*buf != 0 && !isspace(*buf)) {
            if (i >= FIELD_MAX) {
                fprintf(stderr, "Field Name too long \"%s\" \n", buf);
                fclose(input);
                exit(1);
            }
            bed_fld[3][i++] = *buf++;
        }
        bed_fld[3][i] = 0;
        while (isspace(*buf))
            buf++;
        /* Score */
        i = 0;
        while (isdigit(*buf) || *buf == '+' || *buf == '-' || *buf == '.') {
            if (i >= SCORE_MAX) {
                fprintf(stderr, "Score too large \"%s\" \n", buf);
                exit(1);
            }
            bed_fld[4][i++] = *buf++;
        }
        bed_fld[4][i] = 0;
        score = atoi(bed_fld[4]);
        while (isspace(*buf))
            buf++;
        /* Field 6/Strand */
        i = 0;
        while (*buf != 0 && !isspace(*buf)) {
            if (i >= FIELD_MAX) {
                fprintf(stderr, "Field 6 too long \"%s\" \n", buf);
                fclose(input);
                exit(1);
            }
            bed_fld[5][i++] = *buf++;
        }
        bed_fld[5][i] = 0;
        while (isspace(*buf))
            buf++;
        /* Field 7 */
        i = 0;
        while (*buf != 0 && !isspace(*buf)) {
            if (i >= FIELD_MAX) {
                fprintf(stderr, "Field 7 too long \"%s\" \n", buf);
                fclose(input);
                exit(1);
            }
            bed_fld[6][i++] = *buf++;
        }
        bed_fld[6][i] = 0;
        while (isspace(*buf))
            buf++;
        /* Field 8 */
        i = 0;
        while (*buf != 0 && !isspace(*buf)) {
            if (i >= FIELD_MAX) {
                fprintf(stderr, "Field 8 too long \"%s\" \n", buf);
                fclose(input);
                exit(1);
            }
            bed_fld[7][i++] = *buf++;
        }
        bed_fld[7][i] = 0;
        while (isspace(*buf))
            buf++;
        /* Field 9 */
        i = 0;
        while (*buf != 0 && !isspace(*buf)) {
            if (i >= FIELD_MAX) {
                fprintf(stderr, "Field 9 too long \"%s\" \n", buf);
                fclose(input);
                exit(1);
            }
            bed_fld[8][i++] = *buf++;
        }
        bed_fld[8][i] = 0;
        while (isspace(*buf))
            buf++;
        /* Field 10 */
        i = 0;
        while (*buf != 0 && !isspace(*buf)) {
            if (i >= FIELD_MAX) {
                fprintf(stderr, "Field 10 too long \"%s\" \n", buf);
                fclose(input);
                exit(1);
            }
            bed_fld[9][i++] = *buf++;
        }
        bed_fld[9][i] = 0;
        while (isspace(*buf))
            buf++;
        /* Field 11 */
        i = 0;
        while (*buf != 0 && !isspace(*buf)) {
            if (i >= FIELD_MAX) {
                fprintf(stderr, "Field 11 too long \"%s\" \n", buf);
                fclose(input);
                exit(1);
            }
            bed_fld[10][i++] = *buf++;
        }
        bed_fld[10][i] = 0;
        while (isspace(*buf))
            buf++;
        /* Field 12 */
        i = 0;
        while (*buf != 0 && !isspace(*buf)) {
            if (i >= FIELD_MAX) {
                fprintf(stderr, "Field 12 too long \"%s\" \n", buf);
                fclose(input);
                exit(1);
            }
            bed_fld[11][i++] = *buf++;
        }
        bed_fld[11][i] = 0;
        while (isspace(*buf))
            buf++;
#ifdef DEBUG
        printf
            (" [%d]   Chr nb: %s  Start: %lu  End: %lu  Name: %s  Score: %d  Strand: %c  Field #7 %s  Field #8 %s  Field #9 %s  Field #10 %s  Field #11 %s  Field #12 %s\n",
             c++, bed_fld[0], start, end, bed_fld[3], score, bed_fld[5],
             bed_fld[6], bed_fld[7], bed_fld[8], bed_fld[9], bed_fld[10],
             bed_fld[11]);
#endif
        /* Set Strand field  */
        if (bed_fld[5][0] == '+' || bed_fld[5][0] == '-')
            strand = bed_fld[5][0];
        else
            strand = '.';
        /* Set Feature/Name field */
        if (options.feature == NULL) {
            if (bed_fld[3][0] != '\0') {
                options.feature =
                    malloc((strlen(bed_fld[3]) + 1) * sizeof(char));
                strcpy(options.feature, bed_fld[3]);
            }
        }
        /* Print out BED line */
        if (options.oformat == 1) {
            /* Rescale Score */
            new_score = rescale_score(score);
            printf("%s\t%lu\t%lu\t%s%d\t%d\t%c\n", bed_fld[0], start, end,
                   options.feature, count, new_score, strand);
            count++;
        } else if (options.oformat == 2) {
            printf("%s\t%lu\t%lu\t%d\n", bed_fld[0], start, end, score);
        }
    }                           /* End of While */
    if (input != stdin) {
        fclose(input);
    }
    return 0;
}

int main(int argc, char *argv[]) {
#ifdef DEBUG
    mcheck(NULL);
    mtrace();
#endif
    FILE *input;
    int i = 0;
    int usage = 0;
    options.oformat = 1;

    static struct option long_options[] = {
        /* These options may or may not set a flag.
           We distinguish them by their indices. */
        { "debug",      no_argument,       0, 'd' },
        { "help",       no_argument,       0, 'h' },
        { "feature",    required_argument, 0, 'f' },
        { "format",     required_argument, 0, 0   },
        { "oformat",    required_argument, 0, 0   },
        { "minscore",   required_argument, 0, 'm' },
        { "maxscore",   required_argument, 0, 'M' },
        { "name",       required_argument, 0, 0   },
        { "desc",       required_argument, 0, 0   },
        { "color",      required_argument, 0, 0   },
        { "autoscale",  required_argument, 0, 0   },
        { "always0",    required_argument, 0, 0   },
        { "wfunction",  required_argument, 0, 0   },
        { "smoothing",  required_argument, 0, 0   },
        { "visibility", required_argument, 0, 0   },
        { 0,            0,                 0, 0   }
    };

    int option_index = 0;

    while (1) {
        int c =
            getopt_long(argc, argv, "dhf:m:M:o:", long_options,
                        &option_index);
        if (c == -1)
            break;
        switch (c) {
        case 'd':
            options.debug = 1;
            break;
        case 'h':
            options.help = 1;
            break;
        case 'f':
            options.feature = optarg;
            break;
        case 'm':
            minScore = atoi(optarg);
            break;
        case 'M':
            maxScore = atoi(optarg);
            break;
        case 'o':
            fprintf(stderr, "WARNING: option -o is deprecated; please use --format instead.\n");
            options.oformat = atoi(optarg);
            break;
        case 0:
            /* This option is to set the annotation track line */
            if (strcmp(long_options[option_index].name, "name") == 0) {
                options.trackName = optarg;
            }
            if (strcmp(long_options[option_index].name, "desc") == 0) {
                options.trackDesc = optarg;
            }
            if (strcmp(long_options[option_index].name, "color") == 0) {
                options.trackColor = optarg;
            }
            if (strcmp(long_options[option_index].name, "autoscale") == 0) {
                options.autoscale = optarg;
            }
            if (strcmp(long_options[option_index].name, "always0") == 0) {
                options.always0 = optarg;
            }
            if (strcmp(long_options[option_index].name, "wfunction") == 0) {
                options.wfunction = optarg;
            }
            if (strcmp(long_options[option_index].name, "smoothing") == 0) {
                options.smoothing = optarg;
            }
            if (strcmp(long_options[option_index].name, "visibility") == 0) {
                options.visibility = optarg;
            }
            if (strcmp(long_options[option_index].name, "format") == 0) {
                options.oformat = atoi(optarg);
            }
            if (strcmp(long_options[option_index].name, "oformat") == 0) {
                fprintf(stderr, "WARNING: option --oformat is deprecated; please use --format instead.\n");
                options.oformat = atoi(optarg);
            }
            break;
        default:
            printf("?? getopt returned character code 0%o ??\n", c);
            usage = 1;
        }
    }
    if (optind > argc || options.help == 1 || usage || minScore == -1
        || maxScore == -1) {
        fprintf(stderr,
                "Usage: %s [options] [-m <min BED score>] [-M <max BED score>] [<] <BED file|stdin>\n"
                "      - version %s\n" "      where options are:\n"
                "  \t\t -d|--debug                Produce Debug information\n"
                "  \t\t -h|--help                 Show this Help text\n"
                "  \t\t -f|--feature <ft>         Set Feature name <ft>\n"
                "  \t\t --format <1|2>            Set output Format (1:BED[Def] 2:BedGraph)\n"
                "  \t\t -o <1|2>                  Deprecated (please use --format)\n"
                "  \t\t --oformat <1|2>           Deprecated (please use --format)\n"
                "  \t\t --name <name>             Set name for track name field [def. name=SGA-feature]\n"
                "  \t\t --desc <desc>             Set track description field [def. desc=\"ChIP-Seq Custom data\"]\n"
                "  \t\t --color <col>             Define the track color in comma-separated RGB values [def. 100,100,100]\n"
                "  \t\t --autoscale  <on|off>     Data viewing paramenter: set auto-scale to UCSC data view [def=OFF]\n"
                "  \t\t --always0    <on|off>     Data viewing paramenter: always include zero [def=OFF]\n"
                "  \t\t --wfunction  <func>       Data viewing paramenter: windowing function [def=mean+whiskers|maximum|mean|minimum]\n"
                "  \t\t --smoothing  <grade>      Data viewing paramenter: smoothing window [def=OFF[0], <grade>=0,2..16]\n"
                "  \t\t --visibility <mode>       Display mode: [def=dense|full|hide] or [def=1|2|3]\n"
                "\n\tConvert BED file to BED format suitable for displaying ChIP-seq peaks.\n"
                " \tTwo formats are supported: 1) BED Track with Score field; 2) BedGraph.\n"
                " \tFor BED track format, only the graphical parameter 'visibility' can be set.\n\n",
                argv[0], VERSION);
        return 1;
    }
    if (argc > optind) {
        if (!strcmp(argv[optind], "-")) {
            input = stdin;
        } else {
            input = fopen(argv[optind], "r");
            if (NULL == input) {
                fprintf(stderr, "Unable to open '%s': %s(%d)\n",
                        argv[optind], strerror(errno), errno);
                exit(EXIT_FAILURE);
            }
            if (options.debug)
                fprintf(stderr, "Processing file %s\n", argv[optind]);
        }
    } else {
        input = stdin;
    }
    if (options.feature != NULL) {
        char *s = options.feature;
        i = 0;
        while (*s != 0 && !isspace(*s)) {
            if (i >= FT_MAX) {
                fprintf(stderr, "Feature Name too long \"%s\" \n",
                        options.feature);
                return 1;
            }
            s++;
        }
    }
    if (options.trackName == NULL) {
        if (options.feature != NULL) {
            if ((options.trackName =
                 malloc((strlen(options.feature) + 1) * sizeof(char))) ==
                NULL) {
                perror("process_sga: malloc Feature");
                exit(1);
            }
            strcpy(options.trackName, options.feature);
        } else {
            if ((options.trackDesc = malloc(24 * sizeof(char))) == NULL) {
                perror("main: malloc trackDesc");
                exit(1);
            }
            options.trackName = malloc(14 * sizeof(char));      // Allocate memory for trackName
            if (options.trackName == NULL) {
                perror("main: malloc trackName");
                exit(1);
            }
            strcpy(options.trackName, "ChIP-SeqPeaks");
        }
    }
    if (options.trackDesc == NULL) {
        if ((options.trackDesc = malloc(32 * sizeof(char))) == NULL) {
            perror("main: malloc trackDesc");
            exit(1);
        }
        strcpy(options.trackDesc, "ChIP-Seq Custom data");
    }
    if (options.trackColor == NULL) {
        if ((options.trackColor = malloc(12 * sizeof(char))) == NULL) {
            perror("main: malloc trackColor");
            exit(1);
        }
        strcpy(options.trackColor, "0,200,100");
    }
    if (options.autoscale == NULL) {
        if ((options.autoscale = malloc(4 * sizeof(char))) == NULL) {
            perror("main: malloc autoscale");
            exit(1);
        }
        strcpy(options.autoscale, "off");
    }
    if (options.always0 == NULL) {
        if ((options.always0 = malloc(4 * sizeof(char))) == NULL) {
            perror("main: malloc always0");
            exit(1);
        }
        strcpy(options.always0, "off");
    }
    if (options.wfunction == NULL) {
        if ((options.wfunction = malloc(5 * sizeof(char))) == NULL) {
            perror("main: malloc wfunction");
            exit(1);
        }
        strcpy(options.wfunction, "mean");
    }
    if (options.smoothing == NULL) {
        if ((options.smoothing = malloc(4 * sizeof(char))) == NULL) {
            perror("main: malloc smoothing");
            exit(1);
        }
        strcpy(options.smoothing, "off");
    }
    if (options.visibility == NULL) {
        if ((options.visibility = malloc(6 * sizeof(char))) == NULL) {
            perror("main: malloc visibility");
            exit(1);
        }
        if (options.oformat == 1) {
            strcpy(options.visibility, "1");
        } else if (options.oformat == 2) {
            strcpy(options.visibility, "dense");
        }
    }
    if (options.debug) {
        fprintf(stderr, " Arguments:\n");
        if (options.feature != NULL)
            fprintf(stderr, " Feature : %s\n", options.feature);
        fprintf(stderr, " Minimal BED Score: %d\n", minScore);
        fprintf(stderr, " maximal BED score: %d\n", maxScore);
        if (options.oformat == 1) {
            fprintf(stderr, " BED Track Display options:\n");
            fprintf(stderr, " Visibility: %s\n", options.visibility);
            fprintf(stderr, " useScore=1\n");
        } else if (options.oformat == 2) {
            fprintf(stderr, " BedGraph Track Display options:\n");
            fprintf(stderr, " -----------------------------\n");
            fprintf(stderr, " Track name: %s\n", options.trackName);
            fprintf(stderr, " Track description: %s\n", options.trackDesc);
            fprintf(stderr, " Track color: %s\n", options.trackColor);
            fprintf(stderr, " Autoscale: %s\n", options.autoscale);
            fprintf(stderr, " Always0: %s\n", options.always0);
            fprintf(stderr, " Wfunction: %s\n", options.wfunction);
            fprintf(stderr, " Smoothing: %s\n", options.smoothing);
            fprintf(stderr, " Visibility: %s\n", options.visibility);
        }
    }
    if (process_bed(input, argv[optind++]) != 0) {
        return 1;
    }
    return 0;
}
