fit2auc.c

Go to the documentation of this file.

/*****************************************************************************/
#include "tpcclibConfig.h"
/*****************************************************************************/
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
/*****************************************************************************/
#include "tpcextensions.h"
#include "tpcift.h"
#include "tpctac.h"
#include "tpcpar.h"
#include "tpctacmod.h"
/*****************************************************************************/
 
/*****************************************************************************/
static char *info[] = {
  "Calculation of PTAC AUC from specified time to infinity, based on parameters",
  "of decaying sum of exponentials.",
  " ",
  "Usage: @P [Options] parfile starttime [aucfile]",
  " ",
  "Options:",
  " -stdoptions", // List standard options like --help, -v, etc
  " ",
  "See also: fit_dexp, fit2dat, extrapol, paucinf, tacln, parcoll",
  " ",
  "Keywords: AUC, input, pharmacokinetics, clearance, elimination rate",
  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;
  char    parfile[FILENAME_MAX], aucfile[FILENAME_MAX];
  double  starttime=nan("");
 
 
  /*
   *  Get arguments
   */
  if(argc==1) {tpcPrintUsage(argv[0], info, stderr); return(1);}
  parfile[0]=aucfile[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;
    fprintf(stderr, "Error: invalid option '%s'.\n", argv[ai]);
    return(1);
  } else break;
 
  TPCSTATUS status; statusInit(&status);
  statusSet(&status, __func__, __FILE__, __LINE__, TPCERROR_OK);
  status.verbose=verbose-3;
  
  /* 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);}
 
  /* Process other arguments, starting from the first non-option */
  if(ai<argc) strlcpy(parfile, argv[ai++], FILENAME_MAX);
  if(ai<argc) {
    starttime=atofVerified(argv[ai]);
    if(starttime>=0.0) {ai++;}
    else {fprintf(stderr, "Error: invalid start time '%s'.\n", argv[ai]); return(1);}
  }
  if(ai<argc) strlcpy(aucfile, argv[ai++], FILENAME_MAX);
  if(ai<argc) {
    fprintf(stderr, "Error: invalid argument '%s'.\n", argv[ai]);
    return(1);
  }
  /* Did we get all the information that we need? */
  if(!(starttime>=0.0)) {
    fprintf(stderr, "Error: missing command-line argument; use option --help\n");
    return(1);
  }
 
 
  /* In verbose mode print arguments and options */
  if(verbose>1) {
    printf("parfile := %s\n", parfile);
    printf("starttime := %g\n", starttime);
    if(aucfile[0]) printf("aucfile := %s\n", aucfile);
    fflush(stdout);
  }
 
 
  /* 
   *  Read the file
   */
  if(verbose>1) printf("reading %s\n", parfile);
  PAR par; parInit(&par);
  if(parRead(&par, parfile, &status)!=TPCERROR_OK) {
    fprintf(stderr, "Error: %s\n", errorMsg(status.error));
    parFree(&par); return(2);
  }
  if(verbose>2) {
    printf("fileformat := %s\n", parFormattxt(par.format));
    printf("parNr := %d\n", par.parNr);
    printf("tacNr := %d\n", par.tacNr);
  }
 
  /* Try to get data units */
  int cunit=UNIT_UNKNOWN, tunit=UNIT_UNKNOWN;
  {
    int i;
 
    i=iftFindKey(&par.h, "unit", 0);
    if(i<0) i=iftFindKey(&par.h, "calibration_unit", 0);
    if(i>=0) cunit=unitIdentify(par.h.item[i].value);
 
    i=iftFindKey(&par.h, "timeunit", 0);
    if(i<0) i=iftFindKey(&par.h, "time_unit", 0);
    if(i>=0) tunit=unitIdentify(par.h.item[i].value);
    //if(!unitIsTime(tunit)) tunit=UNIT_UNKNOWN;
  }
  if(verbose>2) {
    printf("time_unit := %s\n", unitName(tunit));
    printf("calibration_unit := %s\n", unitName(cunit));
  }
 
 
 
  /*
   *  Allocate place for AUCs
   */
  PAR auc; parInit(&auc);
  if(parAllocate(&auc, 1, par.tacNr)!=TPCERROR_OK) {
    fprintf(stderr, "Error: cannot allocate memory.\n");
    parFree(&par); return(5);
  }
  auc.parNr=1; auc.tacNr=0;
  /* Set the file format */
  auc.format=PAR_FORMAT_UNKNOWN;
  if(aucfile[0]) auc.format=parFormatFromExtension(aucfile);
  if(auc.format==PAR_FORMAT_UNKNOWN) auc.format=PAR_FORMAT_TSV_UK;
  /* Set header contents */
  iftPut(&auc.h, "datafile", parfile, 0, NULL);
  {
    char buf[256];
    time_t t=time(NULL);
    iftPut(&auc.h, "analysis_time", ctime_r_int(&t, buf), 0, NULL);
    tpcProgramName(argv[0], 1, 1, buf, 256);
    iftPut(&auc.h, "program", buf, 0, NULL);
  }
  /* Set parameter name and unit */
  strcpy(auc.n[0].name, "AUC");
  if(unitIsTime(tunit)) {
    int aucunit=unitMultiply(cunit, tunit);
    if(verbose>2) printf("auc_unit := %s\n", unitName(aucunit));
    auc.n[0].unit=aucunit;
  }
 
 
  /*
   *  Compute the AUCs (if possible)
   */
  for(int ti=0; ti<par.tacNr; ti++) {
    if(verbose>2) printf("%s\n", par.r[ti].name);
    double a;
    if(mfEvalIntToInf(modelCode(par.r[ti].model), par.parNr, par.r[ti].p, starttime, &a, verbose-2)) {
      if(verbose>0) fprintf(stderr, "Error: model %s not supported.\n", modelCode(par.r[ti].model)); 
      continue;
    }
    if(verbose>1) printf("AUC for %s := %g\n", par.r[ti].name, a);
    strcpy(auc.r[auc.tacNr].name, par.r[ti].name);
    auc.r[auc.tacNr].p[0]=a;
    auc.tacNr++;
  } // next curve
  if(auc.tacNr==0) {
    fprintf(stderr, "Error: model not supported.\n"); 
    parFree(&par); parFree(&auc);
    return(3);
  }
 
 
  /*
   *  Print and save results
   */
  if(verbose>0 || !aucfile[0]) parWrite(&auc, stdout, PAR_FORMAT_TSV_UK, 0, NULL);
  if(aucfile[0]) {
    /* Save file */
    if(verbose>1) printf("  saving %s\n", aucfile);
    FILE *fp;
    fp=fopen(aucfile, "w");
    if(fp==NULL) {
      fprintf(stderr, "Error: cannot open file for writing.\n");
      parFree(&par); parFree(&auc); return(11);
    }
    int ret=parWrite(&auc, fp, PAR_FORMAT_UNKNOWN, 1, &status);
    fclose(fp);
    if(ret!=TPCERROR_OK) {
      fprintf(stderr, "Error: %s\n", errorMsg(status.error));
      parFree(&par); parFree(&auc); return(12);
    }
    if(verbose>0) printf("AUC(s) saved in %s\n", aucfile);
  }
 
 
  /* Free memory */
  parFree(&par); parFree(&auc);
 
  return(0);
}
/*****************************************************************************/
 
/*****************************************************************************/