taccut.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 "tpctac.h"
/*****************************************************************************/
 
/*****************************************************************************/
static char *info[] = {
  "Extract specified time interval from PET time-activity curves (TACs).",
  " ",
  "Usage: @P [Options] tacfile starttime endtime [outputfile]",
  " ",
  "Options:",
  " -LS | -UptoPeak | -FromPeak | -UptoPeakN | -FromPeakP",
  "     Extract the last sample, samples up to peak, from peak to the end",
  "     up to peak plus one, or from one sample before peak to the end;",
  "     times are not effective with these options.",
  " -stdoptions", // List standard options like --help, -v, etc
  " ",
  "See also: taccat, tactime, tacunit, tacframe, interpol, dftrmovl, imgdelfr",
  " ",
  "Keywords: TAC, tool, simulation, time, cropping",
  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    origNr, ret;
  char   dfile1[FILENAME_MAX], dfile2[FILENAME_MAX];
  TAC    tac, tac2;
  double startT, endT;
  int    calcMode=0; // 0=time range; 1=last sample; 2=up to peak, 3=from peak
 
  /*
   *  Get arguments
   */
  if(argc==1) {tpcPrintUsage(argv[0], info, stderr); return(1);}
  tacInit(&tac); tacInit(&tac2);
  dfile1[0]=dfile2[0]=(char)0;
  startT=endT=nan("");
  /* 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, "LS")==0) {calcMode=1; continue;}
    if(strcasecmp(cptr, "UPTOPEAK")==0) {calcMode=2; continue;}
    if(strcasecmp(cptr, "FROMPEAK")==0) {calcMode=3; continue;}
    if(strcasecmp(cptr, "UPTOPEAKN")==0) {calcMode=4; continue;}
    if(strcasecmp(cptr, "FROMPEAKP")==0) {calcMode=5; 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-1;
  
  /* 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(dfile1, argv[ai++], FILENAME_MAX);}
  if(ai<argc) {
    if(calcMode==0) {
      if(atofCheck(argv[ai], &startT)!=0) {
        fprintf(stderr, "Error: invalid start time '%s'.\n", argv[ai]);
        return(1);
      }
    }
    ai++;
  }
  if(ai<argc) {
    if(calcMode==0) {
      if(atofCheck(argv[ai], &endT)!=0 || endT<=startT) {
        fprintf(stderr, "Error: invalid end time '%s'.\n", argv[ai]);
        return(1);
      }
    }
    ai++;
  }
  if(ai<argc) {strlcpy(dfile2, argv[ai++], FILENAME_MAX);}
  if(ai<argc) {fprintf(stderr, "Error: too many arguments: '%s'.\n", argv[ai]); return(1);}
 
  /* Check that we got what we need */
  if(!dfile1[0]) {tpcPrintUsage(argv[0], info, stderr); return(1);}
  /* Check time range only if it is needed */
  if(calcMode==0 && (isnan(startT) || isnan(endT) || endT<=startT)) {
    fprintf(stderr, "Error: invalid time range.\n"); return(1);}
 
  /* If output file name was not given, then edit the input file */
  if(!dfile2[0]) strcpy(dfile2, dfile1);
 
  /* In verbose mode print arguments and options */
  if(verbose>1) {
    printf("calcMode := %d\n", calcMode);
    printf("dfile1 := %s\n", dfile1);
    printf("dfile2 := %s\n", dfile2);
    if(!isnan(startT)) printf("startT := %g\n", startT);
    if(!isnan(endT)) printf("endT := %g\n", endT);
  }
 
 
  /* 
   *  Read data
   */
  if(verbose>1) printf("reading %s\n", dfile1);
  ret=tacRead(&tac, dfile1, &status);
  if(ret!=TPCERROR_OK) {
    fprintf(stderr, "Error: %s\n", errorMsg(status.error));
    tacFree(&tac); return(2);
  }
  if(verbose>1) {
    printf("fileformat := %s\n", tacFormattxt(tac.format));
    printf("tacNr := %d\n", tac.tacNr);
    printf("sampleNr := %d\n", tac.sampleNr);
  }
  origNr=tac.sampleNr;
 
  /* Sort by sample times */
  if(tacSortByTime(&tac, &status)!=TPCERROR_OK) {
    fprintf(stderr, "Error: %s\n", errorMsg(status.error));
    tacFree(&tac); return(3);
  }
 
  /*
   *  If peak is to be used, then find it, and get time range based on it
   */
  if(calcMode>=2 && calcMode<=5) {
    int maxi; double maxy;
    ret=tacYRange(&tac, -1, NULL, &maxy, NULL, &maxi, NULL, NULL);
    if(ret) {
      fprintf(stderr, "Error: cannot find TAC peak.\n");
      tacFree(&tac); return(4);
    }
    if(verbose>0) {printf("Peak value %g found at %g\n", maxy, tac.x[maxi]); fflush(stdout);}
    if(calcMode==2) {
      startT=-9.9E+30; 
      if(tac.isframe) endT=tac.x2[maxi]; else endT=tac.x[maxi];
    } else if(calcMode==3) {
      if(tac.isframe) startT=tac.x1[maxi]; else startT=tac.x[maxi];
      endT=+9.9E+30; 
    } else if(calcMode==4) {
      startT=-9.9E+30; 
      if(maxi<tac.sampleNr-1) maxi++;
      if(tac.isframe) endT=tac.x2[maxi]; else endT=tac.x[maxi];
    } else if(calcMode==5) {
      if(maxi>0) maxi--;
      if(tac.isframe) startT=tac.x1[maxi]; else startT=tac.x[maxi];
      endT=+9.9E+30; 
    }
  }
 
 
  /*
   *  Extract samples
   */
  if(calcMode==1) {
    /* Extract the last sample */
    ret=tacExtractSamples(&tac, &tac2, tac.sampleNr-1, tac.sampleNr-1);
  } else {
    /* Extract samples inside the user-specified time range */
    ret=tacExtractRange(&tac, &tac2, startT, endT);
  }
  if(ret!=TPCERROR_OK) {
    fprintf(stderr, "Error: invalid range.\n");
    if(verbose>1) printf("error_msg := %s\n", errorMsg(ret));
    tacFree(&tac); tacFree(&tac2); return(4);
  }
  if(verbose>1) {
    printf("original_frameNr := %d\n", origNr);
    printf("final_frameNr := %d\n", tac2.sampleNr);
  }
  tacFree(&tac);
 
  /* If all samples were extracted, and output file name was not given, then do not save output */
  if(origNr==tac2.sampleNr && strcasecmp(dfile1, dfile2)==0) {
    if(verbose>0) printf("nothing done.\n");
    tacFree(&tac2);
    return(0);
  }
 
  /*
   *  Save the extracted TACs
   */
  if(verbose>1) printf("writing %s\n", dfile2);
  FILE *fp; fp=fopen(dfile2, "w");
  if(fp==NULL) {
    fprintf(stderr, "Error: cannot open file for writing (%s)\n", dfile2);
    tacFree(&tac2); return(11);
  }
  ret=tacWrite(&tac2, fp, TAC_FORMAT_UNKNOWN, 1, &status);
  fclose(fp);
  if(ret!=TPCERROR_OK) {
    fprintf(stderr, "Error: %s\n", errorMsg(status.error));
    tacFree(&tac2); return(12);
  }
  if(verbose>0) printf("%d sample(s) were extracted.\n", tac2.sampleNr);
 
  /* Free memory */
  tacFree(&tac2);
 
  return(0);
}
/*****************************************************************************/
 
/*****************************************************************************/