pardiff.c

Go to the documentation of this file.

/*****************************************************************************/
#include "tpcclibConfig.h"
/*****************************************************************************/
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <strings.h>
#include <math.h>
/*****************************************************************************/
#include "tpcextensions.h"
#include "tpcift.h"
#include "tpcpar.h"
/*****************************************************************************/
 
/*****************************************************************************/
static char *info[] = {
  "Calculate difference between parameters in two files as percent change",
  "in file2 from file1 using equation 100%*(result2-result1)/|result1|.",
  "If file name for %-difference is not given, output is printed",
  "on screen (stdout).",
  " ",
  "Usage: @P [options] <file1> <file2> [%-difference file]",
  " ",
  "Options ",
  " -occ[upance]",
  " -inh[ibition]",
  "     Occupance or inhibition percentage is calculated, that is,",
  "     100%*(result1-result2)/result1.",
  " -dif[ference]",
  "     Difference result2-result1 is calculated.",
  " --force",
  "     Result files are allowed to contain different parameter and TAC names",
  "     and units.",
  " -round",
  "     Differences are rounded to the nearest integer percentage.",
  " -stdoptions", // List standard options like --help, -v, etc
  " ",
  "See also: parrenp, parcoll, parai, parformat, parmatch, paradd, parsort",
  " ",
  "Keywords: parameter, modelling, simulation, occupancy, inhibition",
  0};
/*****************************************************************************/
 
/*****************************************************************************/
/* Turn on the globbing of the command line, since it is disabled by default in
   mingw-w64 (_dowildcard=0); in MinGW32 define _CRT_glob instead, if necessary;
   In Unix&Linux wildcard command line processing is enabled by default. */
/*
#undef _CRT_glob
#define _CRT_glob -1
*/
int _dowildcard = -1;
/*****************************************************************************/
 
/*****************************************************************************/
int main(int argc, char **argv)
{
  int ai, help=0, version=0, verbose=1;
  int ri, pi;
  int use_force=0;
  int diff_type=0; // 0=res2/res1-1; 1=1-res2/res1 (occ); 2=1-res2/res1 (inh);
                   // 3=res2-res1 (dif)
  int rounding=0;  // Rounding of percentages: 0=no; 1=yes
  PAR par1, par2;
  char parfile1[FILENAME_MAX], parfile2[FILENAME_MAX], diffile[FILENAME_MAX];
  double f;
 
 
  /*
   *  Get arguments
   */
  if(argc==1) {tpcPrintUsage(argv[0], info, stderr); return(1);}
  parInit(&par1); parInit(&par2);
  parfile1[0]=parfile2[0]=diffile[0]=(char)0;
  /* Options */
  for(ai=1; ai<argc; ai++) if(*argv[ai]=='-') {
    if(tpcProcessStdOptions(argv[ai], &help, &version, &verbose)==0) continue;
    char *cptr=argv[ai]+1; if(*cptr=='-') cptr++; if(!*cptr) continue;
    if(strcasecmp(cptr, "F")==0 || strcasecmp(cptr, "FORCE")==0) {
      use_force=1; continue;
    } else if(strncasecmp(cptr, "NOCHECK", 4)==0) { // deprecated
      use_force=1; continue;
    } else if(strncasecmp(cptr, "OCCUPANCE", 3)==0) {
      diff_type=1; continue;
    } else if(strncasecmp(cptr, "INHIBITION", 3)==0) {
      diff_type=2; continue;
    } else if(strncasecmp(cptr, "DIFFERENCE", 3)==0) {
      diff_type=3; continue;
    } else if(strcasecmp(cptr, "ROUND")==0) {
      rounding=1; continue;
    }
    fprintf(stderr, "Error: invalid option '%s'.\n", argv[ai]); 
    return(1);
  } else break;
 
  /* Print help or version? */
  if(help==2) {tpcHtmlUsage(argv[0], info, ""); return(0);}
  if(help) {tpcPrintUsage(argv[0], info, stdout); return(0);}
  if(version) {tpcPrintBuild(argv[0], stdout); return(0);}
 
  /* Arguments */
  for(; ai<argc; ai++) {
    if(!parfile1[0]) {
      strlcpy(parfile1, argv[ai], FILENAME_MAX); continue;
    } else if(!parfile2[0]) {
      strlcpy(parfile2, argv[ai], FILENAME_MAX); continue;
    } else if(!diffile[0]) {
      strlcpy(diffile, argv[ai], FILENAME_MAX); continue;
    }
    fprintf(stderr, "Error: invalid argument '%s'.\n", argv[ai]);
    return(1);
  }
 
  TPCSTATUS status; statusInit(&status);
  statusSet(&status, __func__, __FILE__, __LINE__, TPCERROR_OK);
  status.verbose=verbose-1;
 
  /* Is something missing? */
  if(!parfile2[0]) {tpcPrintUsage(argv[0], info, stdout); return(1);}
 
  /* In verbose mode print arguments and options */
  if(verbose>1) {
    for(ai=0; ai<argc; ai++) printf("%s ", argv[ai]); 
    printf("\n");
    printf("parfile1 := %s\n", parfile1);
    printf("parfile2 := %s\n", parfile2);
    printf("diffile := %s\n", diffile);
    printf("diff_type := %d\n", diff_type);
    printf("use_force := %d\n", use_force);
    printf("rounding := %d\n", rounding);
    fflush(stdout);
  }
 
 
  /*
   *  Read parameter files
   */
  if(verbose>1) fprintf(stdout, "reading %s\n", parfile1);
  if(parRead(&par1, parfile1, &status)) {
    fprintf(stderr, "Error: %s (%s)\n", errorMsg(status.error), parfile1);
    parFree(&par1); parFree(&par2);
    return(2);
  }
  if(verbose>1) printf("reading %s\n", parfile2);
  if(parRead(&par2, parfile2, &status)) {
    fprintf(stderr, "Error: %s (%s)\n", errorMsg(status.error), parfile2);
    parFree(&par1); parFree(&par2);
    return(3);
  }
 
  /* Sort regions by region name */
  if(verbose>1) fprintf(stdout, "sorting regions\n");
  parSortByName(&par1, &status);
  parSortByName(&par2, &status);
 
  /* Check that contents are matching */
  if(verbose>1) fprintf(stdout, "check for matching results\n");
  if(par1.tacNr!=par2.tacNr) {
    fprintf(stderr, "Error: different TAC number in parameter files.\n");
    parFree(&par1); parFree(&par2); return(4);
  }
  if(par1.parNr!=par2.parNr && use_force==0) {
    fprintf(stderr, "Error: different parameters in parameter files.\n");
    parFree(&par1); parFree(&par2); return(4);
  }
  if(use_force==0 && parCompareTacNames(&par1, &par2, -1, 0, &status)) {
    fprintf(stderr, "Error: different TAC names in parameter files.\n");
    parFree(&par1); parFree(&par2); return(4);
  }
 
  /*
   *  Calculate the differences
   */
  if(verbose>1) fprintf(stdout, "calculating differences\n");
  for(ri=0; ri<par1.tacNr; ri++) {
    for(pi=0; pi<par1.parNr && pi<par2.parNr; pi++) {
      if(par1.r[ri].p[pi]!=0.0) {
        if(diff_type==0) {
          f=par2.r[ri].p[pi]-par1.r[ri].p[pi];
          f/=fabs(par1.r[ri].p[pi]);
          f*=100.0; if(rounding) f=roundf(f);
        } else if(diff_type==3) {
          f=par2.r[ri].p[pi]-par1.r[ri].p[pi];
        } else { // 1 or 2
          f=par1.r[ri].p[pi]-par2.r[ri].p[pi];
          f/=par1.r[ri].p[pi];
          f*=100.0; if(rounding) f=roundf(f);
        }
        par1.r[ri].p[pi]=f;
      } else
        par1.r[ri].p[pi]=0.0;
      /* Set SDs and CLs to NA */
      par1.r[ri].sd[pi]=par1.r[ri].cl1[pi]=par1.r[ri].cl2[pi]=nan("");
    }
  }
  if(par2.parNr<par1.parNr) par1.parNr=par2.parNr;
 
  /*
   *  Write the differences
   */
  if(verbose>1) fprintf(stdout, "writing differences in %s\n", diffile);
  {
    int i, j;
    char buf[256];
    /* change program name */
    tpcProgramName(argv[0], 1, 1, buf, 256);
    i=iftFindKey(&par1.h, "program", 0); if(i>=0) iftDelete(&par1.h, i);
    iftPut(&par1.h, "program", buf, 0, NULL);
    /* change file names */
    i=iftFindKey(&par1.h, "plasmafile", 0); if(i>=0) iftDelete(&par1.h, i);
    i=iftFindKey(&par1.h, "plasmafile2", 0); if(i>=0) iftDelete(&par1.h, i);
    i=iftFindKey(&par1.h, "bloodfile", 0); if(i>=0) iftDelete(&par1.h, i);
    i=iftFindKey(&par1.h, "refname", 0); if(i>=0) iftDelete(&par1.h, i);
    i=iftFindKey(&par1.h, "reffile", 0); if(i>=0) iftDelete(&par1.h, i);
    i=iftFindKey(&par1.h, "datafile", 0); if(i>=0) iftDelete(&par1.h, i);
    iftPut(&par1.h, "reffile", parfile1, 0, NULL);
    iftPut(&par1.h, "datafile", parfile2, 0, NULL);
    /* set study id based on file 2 */
    i=iftFindKey(&par1.h, "studynr", 0);
    if(i<0) i=iftFindKey(&par1.h, "study_number", 0);
    j=iftFindKey(&par2.h, "studynr", 0);
    if(j<0) j=iftFindKey(&par2.h, "study_number", 0);
    if(j>=0) {
      if(i>=0) iftDelete(&par1.h, i);
      iftPut(&par1.h, "studynr", par2.h.item[j].value, 0, NULL);
    }
  }
  /* set units */
  if(diff_type!=3) {
    for(pi=0; pi<par1.parNr; pi++) par1.n[pi].unit=UNIT_PERCENTAGE;
  }
  /* remove WSS */
  for(ri=0; ri<par1.tacNr; ri++) par1.r[ri].wss=nan("");
  /* Save file */
  FILE *fp;
  int ret;
  if(!diffile[0]) fp=stdout; else fp=fopen(diffile, "w");
  if(fp==NULL) {
    fprintf(stderr, "Error: cannot open file for writing.\n");
    parFree(&par1); parFree(&par2); return(11);
  }
  ret=parWrite(&par1, fp, PAR_FORMAT_UNKNOWN, 1, &status);
  if(diffile[0]) fclose(fp); 
  parFree(&par1); parFree(&par2);
  if(ret!=TPCERROR_OK) {
    fprintf(stderr, "Error: %s\n", errorMsg(status.error));
    return(12);
  }
  if(verbose>0 && diffile[0]) 
    printf("parameter differences saved in %s\n", diffile);
  return(0);
}
/*****************************************************************************/
 
/*****************************************************************************/