snaphu_io.c

phase unwrapping algorithm for SAR interferometry

    fprintf(fp,"\n\n");

    /* print an entry for each run-time parameter */
    /* input and output files and main runtime options */
    fprintf(fp,"# File input and output and runtime options\n");
    LogStringParam(fp,"INFILE",infiles->infile);
    fprintf(fp,"LINELENGTH  %ld\n",linelen);
    LogStringParam(fp,"OUTFILE",outfiles->outfile);
    LogStringParam(fp,"WEIGHTFILE",infiles->weightfile);
    if(params->amplitude){
      if(strlen(infiles->ampfile2)){
	LogStringParam(fp,"AMPFILE1",infiles->ampfile);
	LogStringParam(fp,"AMPFILE2",infiles->ampfile2);
      }else{
	LogStringParam(fp,"AMPFILE",infiles->ampfile);
      }
    }else{
      if(strlen(infiles->ampfile2)){
	LogStringParam(fp,"PWRFILE1",infiles->ampfile);
	LogStringParam(fp,"PWRFILE2",infiles->ampfile2);
      }else{
	LogStringParam(fp,"PWRFILE",infiles->ampfile);
      }
    }
    LogStringParam(fp,"MAGFILE",infiles->magfile);
    LogStringParam(fp,"CORRFILE",infiles->corrfile);
    LogStringParam(fp,"ESTIMATEFILE",infiles->estfile);
    LogStringParam(fp,"COSTINFILE",infiles->costinfile);
    LogStringParam(fp,"COSTOUTFILE",outfiles->costoutfile);
    LogStringParam(fp,"LOGFILE",outfiles->logfile);
    if(params->costmode==TOPO){
      fprintf(fp,"STATCOSTMODE  TOPO\n");
    }else if(params->costmode==DEFO){
      fprintf(fp,"STATCOSTMODE  DEFO\n");
    }else if(params->costmode==SMOOTH){
      fprintf(fp,"STATCOSTMODE  SMOOTH\n");
    }else if(params->costmode==NOSTATCOSTS){
      fprintf(fp,"STATCOSTMODE  NOSTATCOSTS\n");
    }
    LogBoolParam(fp,"INITONLY",params->initonly);
    LogBoolParam(fp,"UNWRAPPED_IN",params->unwrapped);
    LogBoolParam(fp,"DEBUG",params->dumpall);
    if(params->initmethod==MSTINIT){
      fprintf(fp,"INITMETHOD  MST\n");
    }else if(params->initmethod==MCFINIT){
      fprintf(fp,"INITMETHOD  MCF\n");
    }
    LogBoolParam(fp,"VERBOSE",params->verbose);

    /* file formats */
    fprintf(fp,"\n# File Formats\n");
    LogFileFormat(fp,"INFILEFORMAT",infiles->infileformat);
    LogFileFormat(fp,"OUTFILEFORMAT",outfiles->outfileformat);
    LogFileFormat(fp,"AMPFILEFORMAT",infiles->ampfileformat);
    LogFileFormat(fp,"MAGFILEFORMAT",infiles->magfileformat);
    LogFileFormat(fp,"CORRFILEFORMAT",infiles->corrfileformat);
    LogFileFormat(fp,"ESTFILEFORMAT",infiles->estfileformat);
    LogFileFormat(fp,"UNWRAPPEDINFILEFORMAT",infiles->unwrappedinfileformat);

    /* SAR and geometry parameters */
    fprintf(fp,"\n# SAR and Geometry Parameters\n");
    fprintf(fp,"ALTITUDE  %.8f\n",
	    params->orbitradius-params->earthradius);
    fprintf(fp,"# ORBITRADIUS  %.8f\n",params->orbitradius);
    fprintf(fp,"EARTHRADIUS  %.8f\n",params->earthradius);
    if(params->bperp){
      fprintf(fp,"BPERP  %.8f\n",params->bperp);
    }else{
      fprintf(fp,"BASELINE %.8f\n",params->baseline);
      fprintf(fp,"BASELINEANGLE_DEG %.8f\n",
	      params->baselineangle*(180.0/PI));
    }
    if(params->transmitmode==PINGPONG){
      fprintf(fp,"TRANSMITMODE  REPEATPASS\n");
    }else if(params->transmitmode==SINGLEANTTRANSMIT){
      fprintf(fp,"TRANSMITMODE  SINGLEANTENNATRANSMIT\n");
    }
    fprintf(fp,"NEARRANGE  %.8f\n",params->nearrange);
    fprintf(fp,"DR  %.8f\n",params->dr);
    fprintf(fp,"DA  %.8f\n",params->da);
    fprintf(fp,"RANGERES  %.8f\n",params->rangeres);
    fprintf(fp,"AZRES  %.8f\n",params->azres);
    fprintf(fp,"LAMBDA  %.8f\n",params->lambda);
    fprintf(fp,"NLOOKSRANGE  %ld\n",params->nlooksrange);
    fprintf(fp,"NLOOKSAZ  %ld\n",params->nlooksaz);
    fprintf(fp,"NLOOKSOTHER  %ld\n",params->nlooksother);
    fprintf(fp,"NCORRLOOKS  %.8f\n",params->ncorrlooks);
    fprintf(fp,"NCORRLOOKSRANGE  %ld\n",params->ncorrlooksrange);
    fprintf(fp,"NCORRLOOKSAZ  %ld\n",params->ncorrlooksaz);
      
    /* scattering model parameters */
    fprintf(fp,"\n# Scattering model parameters\n");
    fprintf(fp,"KDS  %.8f\n",params->kds);
    fprintf(fp,"SPECULAREXP  %.8f\n",params->specularexp);
    fprintf(fp,"DZRCRITFACTOR  %.8f\n",params->dzrcritfactor);
    LogBoolParam(fp,"SHADOW",params->shadow);
    fprintf(fp,"DZEIMIN  %.8f\n",params->dzeimin);
    fprintf(fp,"LAYWIDTH  %ld\n",params->laywidth);
    fprintf(fp,"LAYMINEI  %.8f\n",params->layminei);
    fprintf(fp,"SLOPERATIOFACTOR  %.8f\n",params->sloperatiofactor);
    fprintf(fp,"SIGSQEI  %.8f\n",params->sigsqei);
    
    /* decorrelation model paramters */
    fprintf(fp,"\n# Decorrelation model parameters\n");
    fprintf(fp,"DRHO  %.8f\n",params->drho);
    fprintf(fp,"RHOSCONST1  %.8f\n",params->rhosconst1);
    fprintf(fp,"RHOSCONST2  %.8f\n",params->rhosconst2);
    fprintf(fp,"CSTD1  %.8f\n",params->cstd1);
    fprintf(fp,"CSTD2  %.8f\n",params->cstd2);
    fprintf(fp,"CSTD3  %.8f\n",params->cstd3);
    fprintf(fp,"DEFAULTCORR  %.8f\n",params->defaultcorr);
    fprintf(fp,"RHOMINFACTOR  %.8f\n",params->rhominfactor);
      
    /* PDF model paramters */
    fprintf(fp,"\n# PDF model parameters\n");
    fprintf(fp,"DZLAYPEAK  %.8f\n",params->dzlaypeak);
    fprintf(fp,"AZDZFACTOR  %.8f\n",params->azdzfactor);
    fprintf(fp,"DZEIFACTOR  %.8f\n",params->dzeifactor);
    fprintf(fp,"DZEIWEIGHT  %.8f\n",params->dzeiweight);
    fprintf(fp,"DZLAYFACTOR  %.8f\n",params->dzlayfactor);
    fprintf(fp,"LAYCONST  %.8f\n",params->layconst);
    fprintf(fp,"LAYFALLOFFCONST  %.8f\n",params->layfalloffconst);
    fprintf(fp,"SIGSQSHORTMIN  %ld\n",params->sigsqshortmin);
    fprintf(fp,"SIGSQLAYFACTOR  %.8f\n",params->sigsqlayfactor);

    /* deformation mode paramters */
    fprintf(fp,"\n# Deformation mode parameters\n");
    fprintf(fp,"DEFOAZDZFACTOR  %.8f\n",params->defoazdzfactor);
    fprintf(fp,"DEFOTHRESHFACTOR  %.8f\n",params->defothreshfactor);
    fprintf(fp,"DEFOMAX_CYCLE  %.8f\n",params->defomax);
    fprintf(fp,"SIGSQCORR  %.8f\n",params->sigsqcorr);
    fprintf(fp,"DEFOCONST  %.8f\n",params->defolayconst);

    /* algorithm parameters */
    fprintf(fp,"\n# Algorithm parameters\n");
    fprintf(fp,"INITMAXFLOW  %ld\n",params->initmaxflow);
    fprintf(fp,"ARCMAXFLOWCONST  %ld\n",params->arcmaxflowconst);
    fprintf(fp,"MAXFLOW  %ld\n",params->maxflow);
    fprintf(fp,"KROWEI  %ld\n",params->krowei);
    fprintf(fp,"KCOLEI  %ld\n",params->kcolei);
    fprintf(fp,"KPARDPSI  %ld\n",params->kpardpsi);
    fprintf(fp,"KPERPDPSI  %ld\n",params->kperpdpsi);
    fprintf(fp,"THRESHOLD  %.8f\n",params->threshold);
    fprintf(fp,"INITDZR  %.8f\n",params->initdzr);
    fprintf(fp,"INITDZSTEP  %.8f\n",params->initdzstep);
    fprintf(fp,"MAXCOST  %.8f\n",params->maxcost);
    fprintf(fp,"COSTSCALE  %.8f\n",params->costscale);
    fprintf(fp,"COSTSCALEAMBIGHT  %.8f\n",params->costscaleambight);
    fprintf(fp,"DNOMINCANGLE  %.8f\n",params->dnomincangle);
    fprintf(fp,"NSHORTCYCLE  %ld\n",params->nshortcycle);
    fprintf(fp,"MAXNEWNODECONST  %.8f\n",params->maxnewnodeconst);
    if(params->maxnflowcycles==USEMAXCYCLEFRACTION){
      fprintf(fp,"MAXCYCLEFRACTION  %.8f\n",params->maxcyclefraction);
    }else{
      fprintf(fp,"MAXNFLOWCYCLES  %ld\n",params->maxnflowcycles);
    }
    fprintf(fp,"SOURCEMODE  %ld\n",params->sourcemode);
    fprintf(fp,"CS2SCALEFACTOR  %ld\n",params->cs2scalefactor);
      
    /* file names for dumping intermediate arrays */
    fprintf(fp,"\n# File names for dumping intermediate arrays\n");
    LogStringParam(fp,"INITFILE",outfiles->initfile);
    LogStringParam(fp,"FLOWFILE",outfiles->flowfile);
    LogStringParam(fp,"EIFILE",outfiles->eifile);
    LogStringParam(fp,"ROWCOSTFILE",outfiles->rowcostfile);
    LogStringParam(fp,"COLCOSTFILE",outfiles->colcostfile);
    LogStringParam(fp,"MSTROWCOSTFILE",outfiles->mstrowcostfile);
    LogStringParam(fp,"MSTCOLCOSTFILE",outfiles->mstcolcostfile);
    LogStringParam(fp,"MSTCOSTSFILE",outfiles->mstcostsfile);
    LogStringParam(fp,"RAWCORRDUMPFILE",outfiles->rawcorrdumpfile);
    LogStringParam(fp,"CORRDUMPFILE",outfiles->corrdumpfile);

    /* piece extraction parameters */
    if(params->ntilerow==1 && params->ntilecol==1){
      fprintf(fp,"\n# Piece extraction parameters\n");
      fprintf(fp,"PIECEFIRSTROW  %ld\n",params->piecefirstrow+1);
      fprintf(fp,"PIECEFIRSTCOL  %ld\n",params->piecefirstcol+1);
      fprintf(fp,"PIECENROW  %ld\n",params->piecenrow);
      fprintf(fp,"PIECENCOL  %ld\n",params->piecencol);
    }else{
      fprintf(fp,"\n# Piece extraction parameters\n");
      fprintf(fp,"# Parameters ignored because of tile mode\n");
      fprintf(fp,"# PIECEFIRSTROW  %ld\n",params->piecefirstrow);
      fprintf(fp,"# PIECEFIRSTCOL  %ld\n",params->piecefirstcol);
      fprintf(fp,"# PIECENROW  %ld\n",params->piecenrow);
      fprintf(fp,"# PIECENCOL  %ld\n",params->piecencol);
    }


    /* tile control */
    fprintf(fp,"\n# Tile control\n");
    fprintf(fp,"NTILEROW  %ld\n",params->ntilerow);
    fprintf(fp,"NTILECOL  %ld\n",params->ntilecol);
    fprintf(fp,"ROWOVRLP  %ld\n",params->rowovrlp);
    fprintf(fp,"COLOVRLP  %ld\n",params->colovrlp);
    fprintf(fp,"NPROC  %ld\n",params->nthreads);
    fprintf(fp,"TILECOSTTHRESH  %ld\n",params->tilecostthresh);
    fprintf(fp,"MINREGIONSIZE  %ld\n",params->minregionsize);
    fprintf(fp,"TILEEDGEWEIGHT  %.8f\n",params->tileedgeweight);
    fprintf(fp,"SCNDRYARCFLOWMAX  %ld\n",params->scndryarcflowmax);
    LogBoolParam(fp,"RMTMPTILE",params->rmtmptile);
    if(params->assembleonly){
      LogStringParam(fp,"ASSEMBLEONLY",params->tiledir);
    }else{
      fprintf(fp,"ASSEMBLEONLY  FALSE\n");
    }

    /* connected component control */
    fprintf(fp,"\n# Connected component control\n");
    LogStringParam(fp,"CONNCOMPFILE",outfiles->conncompfile);
    LogBoolParam(fp,"REGROWCONNCOMPS",params->regrowconncomps);
    fprintf(fp,"MINCONNCOMPFRAC  %.8f\n",params->minconncompfrac);
    fprintf(fp,"CONNCOMPTHRESH  %ld\n",params->conncompthresh);
    fprintf(fp,"MAXNCOMPS  %ld\n",params->maxncomps);

    /* close the log file */
    fclose(fp);
  }
}


/* function: LogStringParam()
 * --------------------------
 * Writes a line to the log file stream for the given keyword/value 
 * pair.
 */
void LogStringParam(FILE *fp, char *key, char *value){

  /* see if we were passed a zero length value string */
  if(strlen(value)){
    fprintf(fp,"%s  %s\n",key,value);
    fflush(fp);
  }else{
    fprintf(fp,"# Empty value for parameter %s\n",key);
  }
}


/* LogBoolParam()
 * --------------
 * Writes a line to the log file stream for the given keyword/bool
 * pair.
 */
void LogBoolParam(FILE *fp, char *key, signed char boolvalue){

  if(boolvalue){
    fprintf(fp,"%s  TRUE\n",key);
  }else{
    fprintf(fp,"%s  FALSE\n",key);
  }
}

/* LogFileFormat()
 * ---------------
 * Writes a line to the log file stream for the given keyword/
 * file format pair.
 */
void LogFileFormat(FILE *fp, char *key, signed char fileformat){
  
  if(fileformat==COMPLEX_DATA){
    fprintf(fp,"%s  COMPLEX_DATA\n",key);
  }else if(fileformat==FLOAT_DATA){
    fprintf(fp,"%s  FLOAT_DATA\n",key);
  }else if(fileformat==ALT_LINE_DATA){
    fprintf(fp,"%s  ALT_LINE_DATA\n",key);
  }else if(fileformat==ALT_SAMPLE_DATA){
    fprintf(fp,"%s  ALT_SAMPLE_DATA\n",key);
  }
}


/* function: GetNLines() 
 * ---------------------
 * Gets the number of lines of data in the input file based on the file 
 * size.
 */
long GetNLines(infileT *infiles, long linelen){

  FILE *fp;
  long filesize, datasize;

  /* get size of input file in rows and columns */
  if((fp=fopen(infiles->infile,"r"))==NULL){
    fprintf(sp0,"can't open file %s\n",infiles->infile);
    exit(ABNORMAL_EXIT);
  }
  fseek(fp,0,SEEK_END);
  filesize=ftell(fp);
  fclose(fp);
  if(infiles->infileformat==FLOAT_DATA){
    datasize=sizeof(float);
  }else{
    datasize=2*sizeof(float);
  }
  if(filesize % (datasize*linelen)){
    fprintf(sp0,"extra data in file %s (bad linelength?)\n",
	    infiles->infile);
    exit(ABNORMAL_EXIT);
  }
  return(filesize/(datasize*linelen));               /* implicit floor */

}


/* function: WriteOutputFile()
 * ---------------------------
 * Writes the unwrapped phase to the output file specified, in the
 * format given in the parameter structure.
 */
void WriteOutputFile(float **mag, float **unwrappedphase, char *outfile,
		     outfileT *outfiles, long nrow, long ncol){

  if(outfiles->outfileformat==ALT_LINE_DATA){
    WriteAltLineFile(mag,unwrappedphase,outfile,nrow,ncol);
  }else if(outfiles->outfileformat==ALT_SAMPLE_DATA){
    WriteAltSampFile(mag,unwrappedphase,outfile,nrow,ncol);
  }else if(outfiles->outfileformat==FLOAT_DATA){
    Write2DArray((void **)unwrappedphase,outfile,
		 nrow,ncol,sizeof(float));
  }else{
    fprintf(sp0,"WARNING: Illegal format specified for output file\n");
    fprintf(sp0,"         using default floating-point format\n");
    Write2DArray((void **)unwrappedphase,outfile,
		 nrow,ncol,sizeof(float));
  }
}


/* function: OpenOutputFile()
 * --------------------------
 * Opens a file for writing.  If unable to open the file, tries to 
 * open a file in a dump path.  The name of the opened output file
 * is written into the string realoutfile, for which at least 
 * MAXSTRLEN bytes should already be allocated.
 */
FILE *OpenOutputFile(char *outfile, char *realoutfile){

  char path[MAXSTRLEN], basename[MAXSTRLEN], dumpfile[MAXSTRLEN];
  FILE *fp;

  if((fp=fopen(outfile,"w"))==NULL){

    /* if we can't write to the out file, get the file name from the path */
    /* and dump to the default path */
    ParseFilename(outfile,path,basename);
    StrNCopy(dumpfile,DUMP_PATH,MAXSTRLEN);
    strcat(dumpfile,basename);
    if((fp=fopen(dumpfile,"w"))!=NULL){
      fprintf(sp0,"WARNING: Can't write to file %s.  Dumping to file %s\n",
	     outfile,dumpfile);
      StrNCopy(realoutfile,dumpfile,MAXSTRLEN);
    }else{
      fprintf(sp0,"Unable to write to file %s or dump to file %s\nAbort\n",
	     outfile,dumpfile);
      exit(ABNORMAL_EXIT);
    }
  }else{
    StrNCopy(realoutfile,outfile,MAXSTRLEN);
  }
  return(fp);

}


/* function: WriteAltLineFile()
 * ----------------------------
 * Writes magnitude and phase data from separate arrays to file.
 * Data type is float.  For each line of data, a full line of magnitude data
 * is written, then a full line of phase data.  Dumps the file to a 
 * default directory if the file name/path passed in can