Main Page | Namespace List | Class Hierarchy | Alphabetical List | Class List | File List | Class Members | File Members

unwrap.cc

Go to the documentation of this file.
00001 /*
00002  * Copyright (c) 1999-2003 Bert Kampes
00003  * Copyright (c) 1999-2003 Delft University of Technology, The Netherlands
00004  *
00005  * This file is part of Doris, the Delft o-o radar interferometric software.
00006  *
00007  * Doris program is free software; you can redistribute it and/or modify
00008  * it under the terms of the GNU General Public License as published by
00009  * the Free Software Foundation; either version 2 of the License, or
00010  * (at your option) any later version.
00011  *
00012  * Doris is distributed in the hope that it will be useful,
00013  * but WITHOUT ANY WARRANTY; without even the implied warranty of
00014  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00015  * GNU General Public License for more details.
00016  *
00017  * You should have received a copy of the GNU General Public License
00018  * along with this program; if not, write to the Free Software
00019  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
00020  *
00021  * Publications that contain results produced by the Doris software should
00022  * contain an acknowledgment. (For example: The interferometric processing
00023  * was performed using the freely available Doris software package developed
00024  * by the Delft Institute for Earth-Oriented Space Research (DEOS), Delft
00025  * University of Technology, or include a reference to: Bert Kampes and
00026  * Stefania Usai. \"Doris: The Delft Object-oriented Radar Interferometric
00027  * software.\" In: proceedings 2nd ITC ORS symposium, August 1999. (cdrom)).
00028  *
00029  */
00030 /****************************************************************
00031  * $Source: /users/kampes/DEVELOP/DORIS/doris/src/RCS/unwrap.cc,v $     *
00032  * $Revision: 3.16 $                                            *
00033  * $Date: 2005/04/22 11:30:15 $                                 *
00034  * $Author: kampes $                                            *
00035  *                                                              *
00036  * implementation of unwrapping routines                        *
00037  * -unwrapramon (stanford)                                      *
00038  * -snaphu_unwrap (Curtis Cheng, SU)                            *
00039  ****************************************************************/
00040 
00041 #include "matrixbk.hh"
00042 #include "constants.hh"                 // global constants
00043 #include "productinfo.hh"               // my 'products' class
00044 #include "exceptions.hh"                 // my exceptions class
00045 
00046 #include "unwrap.hh"                    // header file
00047 #include "ioroutines.hh"                // error messages
00048 #include "coregistration.hh"            // distributepoints
00049 #include "readinput.hh"                 // "specified" function
00050 #include "conversion.hh"                // pix2range etc.
00051 #include "orbitbk.hh"                   // getxyz etc.
00052 #include "utilities.hh"                 // getBhBv etc.
00053 
00054 #include <cctype>                       // isspace
00055 #include <cstdio>                       // some compilers, remove function
00056 // #include <cmath>                     // rint function
00057 
00058 
00059 
00060 /****************************************************************
00061  *    unwraptreeframon                                          *
00062  *                                                              *
00063  * unwrap with code howard from ramon:                          *
00064  * this will not work if executable treef_ramon is not in       *
00065  * path.                                                        *
00066  * Place seeds over IFG, and start unwrapping there.  Compose   *
00067  * a map of totally unwrapped regions.  With matlab program     *
00068  * we can check these regions on consistency and adapt them.    *
00069  * Input:                                                       *
00070  *  -                                                           *
00071  * Output:                                                      *
00072  *  -                                                           *
00073  *                                                              *
00074  *    Bert Kampes, 07-Jun-1999                                  *
00075  ****************************************************************/
00076 void unwraptreeframon(
00077         const input_gen     &generalinput,
00078         const input_unwrap  &unwrapinput,
00079         const productinfo   &interferogram)
00080   {
00081   TRACE_FUNCTION("unwraptreeframon (BK 07-Jun-1999)")
00082 
00083 //  char prog[ONE27]="/home/hanssen/projects/cerroprieto/matlab/manuw/treef_ramon ";
00084   char prog[ONE27]              = "treef_ramon ";
00085   char fileouttreeframon[ONE27] = "ramon.uw";
00086   char dummyline[ONE27];
00087 
00088 
00089   const int32 dL = unwrapinput.deltaLseed;
00090   const int32 dP = unwrapinput.deltaPseed;
00091 //  const int32 dL = 100;
00092 //  const int32 dP = dL;
00093 //  char unwrapinput.foregions[ONE27]="REGIONS";
00094 //  char fileoutuint[ONE27]="UWPHASE";
00095 
00096 
00097 // ______ Number of lines/pixels on file (complex interferogram) ______
00098   const int32 fileliness = (interferogram.win.linehi - 
00099                            interferogram.win.linelo + 1) / 
00100                            interferogram.multilookL;
00101   const int32 filepixels = (interferogram.win.pixhi  - 
00102                            interferogram.win.pixlo  + 1) /
00103                            interferogram.multilookP;
00104 
00105 // ______ Allocate matrices ______
00106 // BAND INTERLEAVED !! output of treef_ramon.f (hgt format...)
00107 //  matrix<complr4>     AMP_PHA(fileliness,filepixels);         // read from file ramon.uw
00108   matrix<real4>         AMP(fileliness,filepixels);             // read from file ramon.uw
00109   matrix<real4>         PHA(fileliness,filepixels);             // read from file ramon.uw
00110   matrix<int16>         REGIONS(fileliness,filepixels);         // initialize to 0
00111 
00112 DEBUG.print("only for now like this, maak matrix allocation with setdata??");
00113 //  matrix<real4>       PHASE(fileliness,filepixels,NaN);       // initialize to NaN
00114   matrix<real4>         PHASE(fileliness,filepixels);
00115   PHASE.setdata(NaN);
00116 
00117 
00118 // ______ Make commandstring ______
00119   char basecmdstring[3*ONE27];
00120   ostrstream basecmdstr(basecmdstring,3*ONE27);                 // to convert int to char
00121   basecmdstr.seekp(0);
00122   basecmdstr << prog << interferogram.file << " " 
00123 //           << filepixels << " " << fileliness << " 1 0  ";    // FIRSTLINE = 0 ! (see source)
00124              << filepixels << " " << fileliness << " 0 0  ";    // no ends!!
00125   const int32 posbasecmdstr = basecmdstr.tellp();
00126 
00127 
00128 // ______ Windows for band interleaved file loading, cmp readhgt ______
00129   const uint dummy = 999999;                    // large to force error if not ok
00130   window winamp(1, fileliness, 1, filepixels);
00131   window winpha(1, fileliness, filepixels+1, 2*filepixels);
00132   window offsetbuffer(1,dummy,1,dummy); // dummy not used in readfile, no offset
00133 
00134   register int32 regionnumber = 0;
00135   register int32 line, pixel;
00136   register int32 npoints = 0;
00137   uint seedL,seedP;
00138   matrix<uint> SEEDS;
00139   if (isspace(unwrapinput.seedfile[0]))         // no filename specified
00140     {
00141 // ______ First determine size of matrix the easy way (just counting) ______
00142     for (seedL=4; seedL<=fileliness-4; seedL=seedL+dL)
00143       {
00144       for (seedP=4; seedP<=filepixels-4; seedP=seedP+dP)
00145         {
00146         npoints++;
00147         }
00148       }
00149 
00150 // ______ Then fill the matrix ______
00151     SEEDS.resize(npoints,2);
00152     register int32 cnt = 0;
00153     for (seedL=4; seedL<=fileliness-4; seedL=seedL+dL)
00154       {
00155       for (seedP=4; seedP<=filepixels-4; seedP=seedP+dP)
00156         {
00157         SEEDS(cnt,0) = seedL;
00158         SEEDS(cnt,1) = seedP;
00159         cnt++;
00160         }
00161       }
00162     }
00163 
00164   else                                  // read in file + check
00165     {
00166     INFO.print("Using seedfile for unwrapping seeds.");
00167     if (!existed(unwrapinput.seedfile))
00168       {
00169       PRINT_ERROR("UW_SEED seedfile: does not exist.")
00170       throw(input_error);
00171       }
00172     npoints = filelines(unwrapinput.seedfile);          // number of enters
00173     INFO << "Number of points in seedfile (should be an enter after last one): "
00174          << npoints;
00175     INFO.print();
00176     ifstream ifile(unwrapinput.seedfile, ios::in);      // not binary
00177     SEEDS.resize(npoints,2);
00178 #ifdef __DEBUG
00179     bk_assert(ifile,unwrapinput.seedfile,__FILE__,__LINE__);
00180 #endif
00181 
00182 // ______ Read in data ______
00183     for (line=0; line<npoints; line++)          // read in data
00184       {
00185       ifile >> seedL >> seedP;
00186       SEEDS(line,0) = seedL;
00187       SEEDS(line,1) = seedP;
00188       ifile.getline(dummyline,ONE27,'\n');
00189       }
00190 
00191 // ______ Check input (wrapped interferogram starts at (1,1)) ______
00192     if (max(SEEDS.getcolumn(0)) > fileliness)
00193       {
00194       PRINT_ERROR("max line in seedfile exceeds number of lines on file (complex interferogram).")
00195       throw(input_error);
00196       }
00197     if (max(SEEDS.getcolumn(1)) > filepixels)
00198       {
00199       PRINT_ERROR("max pixel in seedfile exceeds number of pixels on file (complex interferogram).")
00200       throw(input_error);
00201       }
00202     } // seedfile specified
00203 
00204 
00205   for (npoints=0; npoints<SEEDS.lines(); npoints++)
00206     {
00207     seedL = SEEDS(npoints,0);
00208     seedP = SEEDS(npoints,1);
00209 
00210 // ______ Check if region already unwrapped ______
00211     INFO << "Seed (l,p): " << seedL << ", " << seedP;
00212     INFO.print();
00213     if (REGIONS(seedL,seedP) == 0)            // not yet unwrapped
00214       {
00215       regionnumber++;
00216       basecmdstr.seekp(posbasecmdstr);
00217       basecmdstr << seedP << " " << seedL << ends;
00218       INFO << "basecmdstring: " << basecmdstring;
00219       INFO.print();
00220 
00221       remove(fileouttreeframon);                        // try to remove always
00222       system(basecmdstring);
00223 
00224 // ______ Read in output of treef, band interleaved ______
00225 // ______ Use amp for check if unwrapped ok.
00226       readfile (AMP, fileouttreeframon, AMP.lines(), winamp, 
00227       offsetbuffer);
00228       readfile (PHA, fileouttreeframon, PHA.lines(), winpha,
00229       offsetbuffer);
00230 
00231 // ______ Loop over all points, fill REGIONS and PHASE matrix ______
00232       for (line=0; line<AMP.lines(); line++)
00233         {
00234         for (pixel=0; pixel<AMP.pixels(); pixel++)
00235           {
00236           if (AMP(line,pixel) != 0)             // unwrapping ok ?
00237             {
00238             REGIONS(line,pixel) = regionnumber;
00239             PHASE(line,pixel)   = PHA(line,pixel);
00240             }
00241           }
00242         }
00243       }
00244     } // all seeds
00245 
00246   if (remove(fileouttreeframon))
00247     WARNING.print("code 101: could not remove ramon.uw file.");
00248 
00249 
00250 // ______ Remove regions that are not ok unwrapped (too small) ______
00251   const int32 SEVENTY = 70;                             // criterium
00252   const int32 maxregions = regionnumber;
00253   int32 numberofregions = regionnumber;
00254   register int32 numberofpixels  = 0;
00255 
00256   for (regionnumber=1; regionnumber<=maxregions; regionnumber++)
00257     {
00258   
00259 // ______ Count number of pixels ______
00260     for (line=0; line<REGIONS.lines(); line++)
00261       {
00262       if (numberofpixels >= SEVENTY) 
00263         break;                                          // break inner loop
00264       for (pixel=0; pixel<REGIONS.pixels(); pixel++)
00265         {
00266         if (REGIONS(line,pixel) == regionnumber)
00267           {
00268           numberofpixels++;
00269           }
00270         }
00271       }
00272   
00273 // ______ If not ok, delete info ______
00274     if (numberofpixels == 0)
00275       {
00276       WARNING.print("this is impossible..., but still i got this message ones???");
00277       }
00278   
00279     else if (numberofpixels < SEVENTY)                  // too small
00280       {
00281       numberofregions--;
00282       for (line=0; line<REGIONS.lines(); line++)
00283         {
00284         for (pixel=0; pixel<REGIONS.pixels(); pixel++)
00285           {
00286           if (REGIONS(line,pixel) == regionnumber)
00287             {
00288             PHASE(line,pixel) = NaN;                    // delete unwrapped phase
00289             REGIONS(line,pixel) = 0;                    // delete regionnumber
00290             }
00291           }
00292         }
00293       }
00294     }   // loop all regions
00295 
00296 
00297 
00298 // ====== Write result to output files ======
00299   ofstream ofileregions;
00300   openfstream(ofileregions,unwrapinput.foregions,generalinput.overwrit);
00301   bk_assert(ofileregions,unwrapinput.foregions,__FILE__,__LINE__);
00302   ofileregions << REGIONS;
00303   ofileregions.close();
00304 
00305   ofstream ofileuint;
00306   openfstream(ofileuint,unwrapinput.fouint,generalinput.overwrit);
00307   bk_assert(ofileuint,unwrapinput.fouint,__FILE__,__LINE__);
00308   ofileuint << PHASE;
00309   ofileuint.close();
00310 
00311 
00312 // ====== Write info to resultfiles ======
00313   ofstream scratchlogfile("scratchlogunwrap", ios::out | ios::trunc);
00314   scratchlogfile << "\n*******************************************************************"
00315                  << "\n**** Start_unwrap                                             *****"
00316                  << "\nData_output_file: \t\t"
00317                  <<  unwrapinput.fouint
00318                  << "\nData_output_format: \t\t"
00319                  << "real4"
00320                  << "\nData_output_file_regions: \t"
00321                  << unwrapinput.foregions 
00322                  << "\nData_output_format: \t\t"
00323                  << "short int (2B)"
00324                  << "\nProgram for unwrappng: \t\t"
00325                  <<  prog
00326                  << "\nOutput program for unwrapping: \t"
00327                  <<  fileouttreeframon;
00328   if (isspace(unwrapinput.seedfile[0]))         // no seedfile specified
00329     {
00330     scratchlogfile << "\nDelta lines for seed: \t\t"
00331                    <<  dL
00332                    << "\nDelta pixels for seed: \t\t"
00333                    <<  dP;
00334     }
00335   else
00336     { 
00337     scratchlogfile << "\nSeeds where read from file: \t"
00338                    <<  unwrapinput.seedfile;
00339     }
00340   scratchlogfile << "\nNumber of patches used: \t"
00341                  <<  numberofregions
00342                  << "\n* End_unwrap:"
00343                  << "\n*******************************************************************\n";
00344   scratchlogfile.close();
00345 
00346   ofstream scratchresfile("scratchresunwrap", ios::out | ios::trunc);
00347   bk_assert(scratchresfile,"unwrap: scratchresunwrap",__FILE__,__LINE__);
00348   scratchresfile << "\n\n*******************************************************************"
00349                  //<< "\n*_Start_unwrap"
00350                  << "\n*_Start_" << processcontrol[pr_i_unwrap]
00351                  << "\n*******************************************************************"
00352                  << "\nData_output_file:                     \t"
00353                  <<  unwrapinput.fouint
00354                  << "\nData_output_format:                   \t"
00355                  << "real4"
00356                  << "\nData_output_file_regions:             \t"
00357                  << unwrapinput.foregions 
00358                  << "\nData_output_format:                   \t"
00359                  << "short int (2B)"
00360                  << "\nFirst_line (w.r.t. original_master):  \t"
00361                  <<  interferogram.win.linelo
00362                  << "\nLast_line (w.r.t. original_master):   \t"
00363                  <<  interferogram.win.linehi
00364                  << "\nFirst_pixel (w.r.t. original_master): \t"
00365                  <<  interferogram.win.pixlo
00366                  << "\nLast_pixel (w.r.t. original_master):  \t"
00367                  <<  interferogram.win.pixhi
00368                  << "\nMultilookfactor_azimuth_direction:    \t"
00369                  <<  interferogram.multilookL
00370                  << "\nMultilookfactor_range_direction:      \t"
00371                  <<  interferogram.multilookP
00372 
00373                  << "\nProgram for unwrappng:                \t"
00374                  <<  prog
00375                  << "\nOutput program for unwrapping:        \t"
00376                  <<  fileouttreeframon;
00377   if (isspace(unwrapinput.seedfile[0]))         // no seedfile specified
00378     {
00379     scratchresfile
00380                  << "\nDelta lines for seed:                 \t"
00381                  <<  dL
00382                  << "\nDelta pixels for seed:                \t"
00383                  <<  dP;
00384     }
00385   else
00386     { 
00387     scratchresfile << "\nSeeds where read from file:      \t"
00388                    <<  unwrapinput.seedfile;
00389     }
00390   scratchresfile << "\nNumber of patches used:               \t"
00391                  <<  numberofregions
00392                  << "\n*******************************************************************"
00393                  << "\n* End_" << processcontrol[pr_i_unwrap] << "_NORMAL"
00394                  << "\n*******************************************************************\n";
00395   scratchresfile.close();
00396 
00397 // ______ Tidy up ______
00398   PROGRESS.print("finished unwraptreeframon.");
00399   } // END unwraptreeframon
00400 
00401 
00402 /****************************************************************
00403  *    snaphu_unwrap                                             *
00404  *                                                              *
00405  * unwrap with snaphu public domain code, system call           *
00406  * this will not work if executable snaphu is not in path       *
00407  *                                                              *
00408  * Input:                                                       *
00409  *  wrapped ifg                                                 *
00410  * Output:                                                      *
00411  *  unwrapped ifg                                               *
00412  *                                                              *
00413  #%// BK 01-Nov-2002
00414  ****************************************************************/
00415 void snaphu_unwrap(
00416         const input_gen     &generalinput,
00417         const input_unwrap  &unwrapinput,
00418         const productinfo   &interferogram,
00419         const slcimage      &master,
00420         const slcimage      &slave,
00421               orbit         &masterorbit,
00422               orbit         &slaveorbit,
00423         const input_ell     &ellips)
00424   {
00425   TRACE_FUNCTION("snaphu_unwrap (Bert Kampes 01-Nov-2002)")
00426   char prog[ONE27]       = "snaphu ";// run this executable
00427   char configfile[ONE27] = "snaphu.conf";// create this file
00428   INFO << "Configuration file for snaphu: " << configfile;
00429   INFO.print();
00430 
00431 
00432   // ______ Number of lines/pixels on file (complex interferogram) ______
00433   const int32 filelines  = (interferogram.win.linehi - 
00434                            interferogram.win.linelo + 1) / 
00435                            interferogram.multilookL;
00436   const int32 filepixels = (interferogram.win.pixhi  - 
00437                            interferogram.win.pixlo  + 1) /
00438                            interferogram.multilookP;
00439 
00440   // ______ Make commandstring ______
00441   char basecmdstring[3*ONE27];
00442   ostrstream basecmdstr(basecmdstring,3*ONE27);// to convert int to char
00443   basecmdstr.seekp(0);
00444   basecmdstr << prog << " -f " << configfile << " "
00445              << interferogram.file << " " 
00446              << filepixels << ends;
00447 
00448 
00449   // _____ Compute some parameters for snaphu here ______
00450   real8 NEARRANGE = master.pix2range(interferogram.win.pixlo);
00451     //pix2range(interferogram.win.pixlo,master.t_range1,master.rsr2x);
00452   real8 DR        = real8(interferogram.multilookP) *
00453     (master.pix2range(interferogram.win.pixlo+1.0) - NEARRANGE);
00454     //(pix2range(interferogram.win.pixlo+1,master.t_range1,master.rsr2x) - NEARRANGE);
00455   real8 RANGERES  = ((master.rsr2x/2.)/master.rbw) * 
00456                     (DR/real8(interferogram.multilookP));
00457 
00458   cn P;   // point at ellips for mid image, returned by lp2xyz
00459   real8 line  = 0.5*real8(interferogram.win.linelo+interferogram.win.linehi);
00460   real8 pixel = 0.5*real8(interferogram.win.pixlo+interferogram.win.pixhi);
00461   lp2xyz(line,pixel,ellips,master,masterorbit,P);// returns P
00462   real8 EARTHRADIUS = P.norm();
00463   // ______ Compute xyz for master satellite ______
00464   real8 time_azi;                               // returned
00465   real8 time_range;                             // returned
00466   xyz2t(time_azi,time_range,master, masterorbit, P);
00467   cn M = masterorbit.getxyz(time_azi);          // master satellite position
00468   xyz2t(time_azi,time_range,slave, slaveorbit, P);
00469   cn S = slaveorbit.getxyz(time_azi);           // slave satellite position
00470   real8 ORBITRADIUS = M.norm();
00471   // ______ Baseline parametrizations returned ______
00472   real8 Bh, Bv, Bpar, Bperp, theta;
00473   real8 BASELINE;
00474   real8 BASELINEANGLE_RAD;
00475   BalphaBhBvBparBperpTheta(
00476     BASELINE, BASELINEANGLE_RAD,
00477     Bh, Bv, Bpar, Bperp, theta, M, P, S);
00478   real8 BASELINEANGLE_DEG = rad2deg(BASELINEANGLE_RAD);
00479   ;// dA via orbit...
00480   lp2xyz(line+1.,pixel,ellips,master,masterorbit,P);// returns P
00481   xyz2t(time_azi,time_range,master, masterorbit, P);// returns time
00482   cn M1 = masterorbit.getxyz(time_azi);
00483   real8 DA        = abs(real8(interferogram.multilookL)*M.dist(M1));
00484   real8 AZRES     = ((master.prf)/master.abw) * 
00485                     (DA/real8(interferogram.multilookL));
00486 
00487   // ______ Create config file ______
00488   ofstream snaphuconfigfile(configfile, ios::out | ios::trunc);
00489   snaphuconfigfile
00490     << "# snaphu configuration file\n"
00491     << "#\n"
00492     << "# Lines with fewer than two fields and lines whose first non-whitespace\n"
00493     << "# characters are not alphnumeric are ignored.  For the remaining lines,\n"
00494     << "# anything after the first two fields (delimited by whitespace) is\n"
00495     << "# also ignored.  Inputs are converted in the order they appear in the file;\n"
00496     << "# if multiple assignments are made to the same parameter, the last one\n"
00497     << "# given is the one used.  Parameters in this file will be superseded by\n"
00498     << "# parameters given on the command line after the -f flag specifying this\n"
00499     << "# file.  Multiple configuration files may be given on the command line.\n"
00500     << "\n"
00501     << "# CONFIG FOR SNAPHU\n"
00502     << "# -----------------------------------------------------------\n"
00503     << "# Created by Doris software#\n"
00504     << "# Has been run from within Doris with command:\n\n"
00505     << "#    " << basecmdstring << "\n"
00506     << "# -----------------------------------------------------------\n"
00507     << "\n"
00508     << "# Statistical-cost mode (TOPO, DEFO, SMOOTH, or NOSTATCOSTS)\n"
00509     << "STATCOSTMODE    " << unwrapinput.snaphu_mode << "\n"
00510     << "\n"
00511     << "# Output file name\n"
00512     << "OUTFILE         " << unwrapinput.fouint << "\n"
00513     << "\n"
00514     << "# Correlation file name\n";
00515   if (specified(unwrapinput.snaphu_coh))
00516     {
00517     DEBUG.print("Using coherence of Doris for snaphu:");
00518     DEBUG.print("Please make sure window sizes used in config file are correct");
00519     DEBUG.print("Maybe you will have to re-run snaphu by hand? (see notes in snaphu.conf)");
00520     snaphuconfigfile
00521       << "CORRFILE      " << unwrapinput.snaphu_coh << "\n"
00522       << "\n"
00523       << "# Correlation file format\n"
00524       << "CORRFILEFORMAT        FLOAT_DATA\n";
00525     }
00526   else
00527     { 
00528     snaphuconfigfile
00529       << "# CORRFILE      " << "doris.coh" << "\n"
00530       << "\n"
00531       << "# Correlation file format\n"
00532       << "# CORRFILEFORMAT        FLOAT_DATA\n";
00533     }
00534   snaphuconfigfile
00535     << "# Text file to which runtime parameters will be logged.  The format of\n"
00536     << "# that file will be suitable so that it can also be used as a\n"
00537     << "# configuration file.\n";
00538   if (specified(unwrapinput.snaphu_log))
00539     {
00540     snaphuconfigfile
00541       << "LOGFILE           " << unwrapinput.snaphu_log
00542       << "\n";
00543     }
00544   else
00545     {
00546     snaphuconfigfile
00547       << "# LOGFILE           " << "snaphu.log"
00548       << "\n";
00549     }
00550   snaphuconfigfile
00551     << "# Algorithm used for initialization of wrapped phase values.  Possible\n"
00552     << "# values are MST and MCF.\n"
00553     << "INITMETHOD      " << unwrapinput.snaphu_init << "\n"
00554     << "\n"
00555     << "# Verbose-output mode (TRUE or FALSE)\n"
00556     << "VERBOSE         " << unwrapinput.snaphu_verbose << "\n"
00557     << "\n"
00558     << "\n"
00559     << "################\n"
00560     << "# File formats #\n"
00561     << "################\n"
00562     << "\n"
00563     << "# Valid data formats:\n"
00564     << "#\n"
00565     << "# COMPLEX_DATA:      complex values: real, imag, real, imag\n"
00566     << "# ALT_LINE_DATA:     real values from different arrays, alternating by line\n"
00567     << "# ALT_SAMPLE_DATA:   real values from different arrays, alternating by sample\n"
00568     << "# FLOAT_DATA:        single array of floating-point data\n"
00569     << "#\n"
00570     << "\n"
00571     << "# Input file format\n"
00572     << "INFILEFORMAT          COMPLEX_DATA\n"
00573     << "\n"
00574     << "# Output file format (this is almost hgt, BK)\n";
00575   if (unwrapinput.oformatflag == FORMATHGT)
00576     {
00577     snaphuconfigfile
00578       << "OUTFILEFORMAT         ALT_LINE_DATA\n"
00579       << "# OUTFILEFORMAT         FLOAT_DATA\n";
00580     }
00581   else if (unwrapinput.oformatflag == FORMATR4)
00582     {
00583     snaphuconfigfile
00584       << "# OUTFILEFORMAT         ALT_LINE_DATA\n"
00585       << "OUTFILEFORMAT         FLOAT_DATA\n";
00586     }
00587   else
00588     {
00589     PRINT_ERROR("unknown format specified for output with SNAPHU")
00590     throw(unhandled_case_error);
00591     }
00592   snaphuconfigfile
00593     << "\n"
00594     << "\n"
00595     << "###############################\n"
00596     << "# SAR and geometry parameters #\n"
00597     << "###############################\n"
00598     << "\n"
00599     << "# Orbital radius (double, meters) or altitude (double, meters).  The\n"
00600     << "# radius should be the local radius if the orbit is not circular.  The\n"
00601     << "# altitude is just defined as the orbit radius minus the earth radius.\n"
00602     << "# Only one of these two parameters should be given.\n"
00603     << "#ORBITRADIUS             7153000.0 (example)\n"
00604     << "ORBITRADIUS             "
00605     << setiosflags(ios::fixed) << setw(10) << setprecision(2)
00606     << ORBITRADIUS << "\n"
00607     << "#ALTITUDE               775000.0\n"
00608     << "\n"
00609     << "# Local earth radius (double, meters).  A spherical-earth model is\n"
00610     << "# used.\n"
00611     << "#EARTHRADIUS             6378000.0 (example)\n"
00612     << "EARTHRADIUS             " 
00613     << setiosflags(ios::fixed) << setw(10) << setprecision(2)
00614     << EARTHRADIUS << "\n"
00615     << "\n"
00616     << "# The baseline parameters are not used in deformation mode, but they\n"
00617     << "# are very important in topography mode.  The parameter BASELINE\n"
00618     << "# (double, meters) is the physical distance (always positive) between\n"
00619     << "# the antenna phase centers.  The along-track componenet of the\n"
00620     << "# baseline is assumed to be zero.  The parameter BASELINEANGLE_DEG\n"
00621     << "# (double, degrees) is the angle between the antenna phase centers\n"
00622     << "# with respect to the local horizontal.  Suppose the interferogram is\n"
00623     << "# s1*conj(s2).  The baseline angle is defined as the angle of antenna2\n"
00624     << "# above the horizontal line extending from antenna1 towards the side\n"
00625     << "# of the SAR look direction.  Thus, if the baseline angle minus the\n"
00626     << "# look angle is less than -pi/2 or greater than pi/2, the topographic\n"
00627     << "# height increases with increasing elevation.  The units of\n"
00628     << "# BASELINEANGLE_RAD are radians.\n"
00629     << "#BASELINE                150.0 (example)\n"
00630     << "#BASELINEANGLE_DEG       225.0 (example)\n"
00631     << "#BASELINEANGLE_RAD      3.92699 (example)\n"
00632     << "BASELINE                " << BASELINE << "\n"
00633     << "BASELINEANGLE_DEG       " << BASELINEANGLE_DEG << "\n"
00634     << "#BASELINEANGLE_RAD      " << BASELINEANGLE_RAD << "\n"
00635     << "\n"
00636     << "# If the BPERP parameter is given, the baseline angle is taken to be\n"
00637     << "# equal to the look angle (mod pi) at midswath, and the length of the\n"
00638     << "# baseline is set accordingly.  Particular attention must be paid to\n"
00639     << "# the sign of this parameter--it should be negative if increasing\n"
00640     << "# phase implies increasing topographic height.\n"
00641     << "#BPERP          -150.0 (example)\n"
00642     << "#BPERP          " << Bperp << "\n"
00643     << "\n"
00644     << "# The transmit mode should be either REPEATPASS or PINGPONG if both\n"
00645     << "# antennas transmitted and both received (REPEATPASS and PINGPONG have\n"
00646     << "# the same effect); the transmit mode should be SINGLEANTENNATRANSMIT\n"
00647     << "# if only one antenna was used to transmit while both antennas\n"
00648     << "# received.  In single-antenna-transmit mode, the baseline is\n"
00649     << "# effectively halved.  This parameter is ignored for cost modes other\n"
00650     << "# than topography.\n"
00651     << "TRANSMITMODE    REPEATPASS\n"
00652     << "\n"
00653     << "# Slant range from platform to first range bin in input data file\n"
00654     << "# (double, meters).  Be sure to modify this parameter if the input\n"
00655     << "# file is extracted from a larger scene.  The parameter does not need\n"
00656     << "# to be modified is snaphu is unwrapping only a subset of the input file.\n"
00657     << "#NEARRANGE       831000.0 (example)\n"
00658     << "NEARRANGE       " 
00659     << setiosflags(ios::fixed) << setw(20) << setprecision(2)
00660     << NEARRANGE << "\n"
00661     << setw(20) << setprecision(10)
00662     << "\n"
00663     << "# Slant range and azimuth pixel spacings of input interferogram after\n"
00664     << "# any multilook averaging.  This is not the same as the resolution.\n"
00665     << "# (double, meters).\n"
00666     << "#DR              8.0 (example)\n"
00667     << "#DA              20.0 (example)\n"
00668     << "DR              " << DR << "\n"
00669     << "DA              " << DA << "\n"
00670     << "\n"
00671     << "# Single-look slant range and azimuth resolutions.  This is not the\n"
00672     << "# same as the pixel spacing.  (double, meters).\n"
00673     << "#RANGERES        10.0 (example)\n"
00674     << "#AZRES           6.0 (example)\n"
00675     << "RANGERES        " << RANGERES << "\n"
00676     << "AZRES           " << AZRES    << "\n"
00677     << "\n"
00678     << "# Wavelength (double, meters).\n"
00679     << "#LAMBDA          0.0565647 (example)\n"
00680     << "LAMBDA          " << master.wavelength << "\n"
00681     << "\n"
00682     << "# Number of real (not necessarily independent) looks taken in range and\n"
00683     << "# azimuth to form the input interferogram (long).\n"
00684     << "NLOOKSRANGE     " << interferogram.multilookP << "\n"
00685     << "NLOOKSAZ        " << interferogram.multilookL << "\n"
00686     << "\n"
00687     << "\n"
00688     << "# Equivalent number of independent looks (double, dimensionless) that were\n"
00689     << "# used to generate correlation file if one is specified.  This parameter\n"
00690     << "# is ignored if the correlation data are generated by the interferogram\n"
00691     << "# and amplitude data.\n"
00692     << "#"
00693     << "# The equivalent number of independent looks is approximately equal to the\n"
00694     << "# real number of looks divided by the product of range and azimuth\n"
00695     << "# resolutions, and multiplied by the product of the single-look range and\n"
00696     << "# azimuth pixel spacings.  It is about 0.53 times the number of real looks\n"
00697     << "# for ERS data processed without windowing.\n"
00698     << "# (BK: taken from Curtis example config file.)\n"
00699     << "NCORRLOOKS      23.8\n"
00700     << "\n"
00701     << "# Number of looks that should be taken in range and azimuth for estimating\n"
00702     << "# the correlation coefficient from the interferogram and the amplitude\n"
00703     << "# data.  These numbers must be larger than NLOOKSRANGE and NLOOKSAZ.\n"
00704     << "# The actual numbers used may be different since we prefer odd integer\n"
00705     << "# multiples of NLOOKSRANGE and NLOOKSAZ (long).  These numbers are ignored\n"
00706     << "# if a separate correlation file is given as input.\n"
00707     << "# (BK: taken from Curtis example config file.)\n"
00708     << "NCORRLOOKSRANGE 3\n"
00709     << "NCORRLOOKSAZ    15\n"
00710     << "\n"
00711     << "# End of snaphu configuration file\n"
00712     << "\n";
00713   snaphuconfigfile.close();
00714   PROGRESS << "Configuration file for snaphu: " << configfile << " created." << ends;
00715   PROGRESS.print();
00716 
00717 
00718 
00719   // ______ Call snaphu! ______
00720   INFO.print("System call for running snaphu follows:");
00721   INFO.print(basecmdstring);
00722   system(basecmdstring);
00723 
00724   // ====== Write info to resultfiles ======
00725   ofstream scratchlogfile("scratchlogunwrap", ios::out | ios::trunc);
00726   scratchlogfile
00727     << "\n*********************************************************"
00728     << "\n**** Start_unwrap                                   *****"
00729     << "\n*********************************************************"
00730     << "\nData_output_file: \t\t"
00731     <<  unwrapinput.fouint
00732     << "\nData_output_format: \t\t"
00733     << "hgt"
00734     << "\nProgram for unwrappng: \t\t"
00735     <<  prog
00736     << "\n* End_unwrap:"
00737     << "\n*********************************************************\n";
00738   scratchlogfile.close();
00739 
00740   ofstream scratchresfile("scratchresunwrap", ios::out | ios::trunc);
00741   bk_assert(scratchresfile,"unwrap: scratchresunwrap",__FILE__,__LINE__);
00742   scratchresfile 
00743     << "\n\n*******************************************************************"
00744     << "\n*_Start_" << processcontrol[pr_i_unwrap]
00745     << "\n*******************************************************************"
00746     << "\nData_output_file:                     \t"
00747     <<  unwrapinput.fouint
00748     << "\nData_output_format:                   \t"
00749     << "hgt"
00750     << "\nFirst_line (w.r.t. original_master):  \t"
00751     <<  interferogram.win.linelo
00752     << "\nLast_line (w.r.t. original_master):   \t"
00753     <<  interferogram.win.linehi
00754     << "\nFirst_pixel (w.r.t. original_master): \t"
00755     <<  interferogram.win.pixlo
00756     << "\nLast_pixel (w.r.t. original_master):  \t"
00757     <<  interferogram.win.pixhi
00758     << "\nMultilookfactor_azimuth_direction:    \t"
00759     <<  interferogram.multilookL
00760     << "\nMultilookfactor_range_direction:      \t"
00761     <<  interferogram.multilookP
00762     << "\nNumber of lines (multilooked):        \t"
00763     <<  filelines 
00764     << "\nNumber of pixels (multilooked):       \t"
00765     <<  filepixels
00766     << "\nProgram used for unwrapping:          \t"
00767     <<  prog
00768     << "\n*******************************************************************"
00769     << "\n* End_" << processcontrol[pr_i_unwrap] << "_NORMAL"
00770     << "\n*******************************************************************\n";
00771   scratchresfile.close();
00772 
00773   // ______ Tidy up ______
00774   PROGRESS.print("Finished snaphu_unwrap.");
00775   } // END snaphu_unwrap
00776 

Generated on Fri Apr 22 15:58:03 2005 for Doris by doxygen 1.3.6