/* button.c

   Compile with:   gcc -g button.c -lm -o buttonx
                   gcc -g button.c -DGAP -lm -o buttong 

   Example run:    button beta_reference S=6

   Initiates agree test runs for a fixed specified min_length (S) and
   values of the minimum required percent identity between 10 and 100,
   in increments of 1%. Both the strict (no gaps allowed, X) and
   general (gaps allowed, G) cases are considered. For each such run 
   it computes the # of false negatives/positives with respect to a
   reference set of landmarks.

   NOTE: Names of alignment files/paths may need to be changed. Look
   for LLL to locate the references.

*/

#include <stdio.h>
#include <stdlib.h>
#include <strings.h>
#include <ctype.h>
#include <sys/types.h>
#include <sys/param.h>


#define Usage "%s <ref_file> [S=<min_region_size>] [L=<lower>] [U=<upper>]\n"

void fatal(char *message);
void fatalf(char *msg, char *val);
FILE *ckopen(char *name, char *mode);
void *ckalloc(size_t amount);

typedef struct data {
       char filename[50];
       int  nfp;
       int  nfn;
} data_t;

void run_agree(int low, int up, int epsilon, int flag,
                 data_t *low_data, data_t *up_data);
void run_test(int percent, data_t *data);
void compare(char *filename, data_t *data);


#define    DEFAULT_LOWER_BOUND    10
#define    DEFAULT_UPPER_BOUND    100
#define    DEFAULT_REFINE         1
#define    DEFAULT_LENGTH         6

#define    LOWER	0x1
#define    UPPER        0x2


int *Fct, *Regs, To, From, Length, len;

int main(int argc, char *argv[]) 
{
    int    to, from, i;
    int    low, up, eps;
    data_t low_data, up_data; 
    char   buffer[200];
    FILE  *fp;

    if ((argc<3) || (argc>5)) fatalf(Usage,argv[0]);

    Length = DEFAULT_LENGTH;
    low = DEFAULT_LOWER_BOUND;
    up = DEFAULT_UPPER_BOUND;
    eps = DEFAULT_REFINE;

    while (argc>2) {
       if (!strncmp(argv[argc-1],"S=",2))
           Length = atoi(argv[argc-1]+2);
       else if (!strncmp(argv[argc-1],"L=",2))
           low = atof(argv[argc-1]+2);
       else if (!strncmp(argv[argc-1],"U=",2))
           up = atof(argv[argc-1]+2);
       else 
           fatalf(Usage,argv[0]);           
       argc--;
    }

    /* read reference file */
    fp = ckopen(argv[1],"r");
    len = 0;
    fgets(buffer,200,fp);
    if (buffer[0]!='>') fatal("Wrong file configuration.");
    sscanf(buffer+1, "%d %d", &From, &To);

 
    Fct = (int *)ckalloc((To-From+1)*sizeof(int))-From;
    Regs = (int *)ckalloc((To-From+1)*sizeof(int))-From;
    for (i=From; i<=To; Fct[i]=Regs[i]=0, i++);

    while (fgets(buffer,200,fp)!=NULL) {
       if (!isdigit(buffer[0])) continue;
       sscanf(buffer, "%d %d", &from, &to);
       for (i=from; i<=to; i++) 
            if ((i>=From) && (i<=To)) Fct[i] = 1;
       len += to-from+1;
    }
    fclose(fp);


    run_agree(low, up, eps, LOWER | UPPER, &low_data, &up_data); 
/*
    system("sort -nr +7 DEVIATION > araCe_DEVIATION.total");
*/
    return 0;
}

void run_agree(int low, int up, int epsilon, int flag, data_t *low_data, data_t *up_data)
{
     int mid;
     data_t  mid_data;
     FILE   *fp;

     if (up-low < epsilon) return;
     if (flag & LOWER) 
         run_test(low,low_data); 
     if (flag & UPPER) 
         run_test(up, up_data);
       
     if ((low_data->nfp!=up_data->nfp) || 
         (low_data->nfn!=up_data->nfn)) {
          fp = ckopen("DEVIATION","a");
          if (flag & LOWER) 
              fprintf(fp, "%-15s   \tfp: %d \tfn: %d   \ttotal: %d (%d)\n",
                       low_data->filename, low_data->nfp, low_data->nfn,
                       low_data->nfp+low_data->nfn, len); 
          fclose(fp);     
          mid = (up+low)>>1; 
          if (mid!=up) 
              run_agree(low,mid,epsilon,UPPER,low_data,&mid_data);
          if (mid!=low)
              run_agree(mid,up,epsilon,LOWER,&mid_data,up_data); 
          if ((flag & UPPER) && ((up_data->nfp!=mid_data.nfp) || 
                                 (up_data->nfn!=mid_data.nfn)))
              fprintf(fp, "%-15s   \tfp: %d \tfn: %d   \ttotal: %d (%d)\n",
                       up_data->filename, up_data->nfp, up_data->nfn,
                       up_data->nfp+up_data->nfn, len);
     } else {
          if (flag & LOWER) { 
              fp = ckopen("DEVIATION","a");
              fprintf(fp, "%-15s   \tfp: %d \tfn: %d   \ttotal: %d (%d)\n",
                       low_data->filename, low_data->nfp, low_data->nfn,
                       low_data->nfp+low_data->nfn, len);
              fclose(fp);
          }
     }
     if ((flag & LOWER) && (flag & UPPER)) {
          fp = ckopen("DEVIATION","a");
          fprintf(fp, "%-15s   \tfp: %d \tfn: %d   \ttotal: %d (%d)\n",
                       up_data->filename, up_data->nfp, up_data->nfn,
                       up_data->nfp+up_data->nfn, len);
          fclose(fp);
     }
}

void run_test(int percent, data_t *data)
{
     char command[200]; 
     char filename[50];

     /* LLL change name of region here */
     sprintf(filename, "agree.araCe.l%d.p%d.list", Length, percent);

     /* LLL - example of command line for tests on the beta-globin alignment 
     sprintf(command,"agree alignment -I /usr/local/globin/data/align.beta -r %d %d -l %d -a 3 -p %d -g > %s", 
                     From, To, Length, percent, filename);
     */

     /* LLL change name of alignment file + path here */
#ifndef GAP
     sprintf(command,"/net/finch/usr/local/globin/bin/agree araCe.a -I /home/galapagos/florea/test.d -r %d %d -l %d -a 5 -p %d -x > %s", 
                     From, To, Length, percent, filename);
#else 
     sprintf(command,"/net/finch/usr/local/globin/bin/agree araCe.a -I /home/galapagos/florea/test.d -r %d %d -l %d -a 5 -p %d -g > %s",  
                     From, To, Length, percent, filename);
#endif     

     system(command); 
     compare(filename, data);

     return;
}

void compare(char *filename, data_t *data)
{
     FILE *tfp;
     int   from, to, i;
     char  buffer[200];

     tfp = ckopen(filename,"r");
     while (fgets(buffer,200,tfp)!=NULL) {
        if (!isdigit(buffer[0])) continue;
        sscanf(buffer, "%d %d", &from, &to);
        for (i=from; i<=to; i++)
             if ((i>=From) && (i<=To)) Regs[i] = 1;
     }       
     fclose(tfp);             

     strcpy(data->filename, filename);
     data->nfn = data->nfp = 0;
     for (i=From; i<=To; i++)
          if (!Fct[i] && Regs[i]) data->nfp++;
          else if (Fct[i] && !Regs[i]) data->nfn++;
     for (i=From; i<=To; Regs[i++]=0);

     return;
}


void fatal(char *msg)
{
        (void)fprintf(stderr, "%s\n", msg);
        exit(1);
}
 
void fatalf(char *msg, char *val)
{
        (void)fprintf(stderr, msg, val);
        (void)putc('\n', stderr);
        exit(1);
}    

FILE *ckopen(char *name, char *mode)
{
        FILE *fp;
 
        if ((fp = fopen(name, mode)) == NULL)
                fatalf("Cannot open %s.", name);
        return fp;
}
        
 
void *ckalloc(size_t amount)
{
        void *p;
        
        if ((p = malloc(amount)) == NULL)
                fatal("Ran out of memory.");
        return p;
}