Home

Delft Object-oriented Radar Interferometric Software User's manual

image

Contents

1. Figure 32 1 Geometric configuration for slant to height conversion 109 32 2 Output Description The process control flag is switched at successful exit slant2height 1l An example of the output in the products result file KKKKKKKKKKKKKKKKKKKKK KK KK ckckckckckckckckckckckckckck ckckckckckckckckckckckckck ck ck ck ck ck ck ck ckckckck k ck kk Start slantzh ck ck Ck ck ckckckckckckckckckckckckckck ck ckckckckckckckckckckckckckckckckckckckckckckckckckckckckckckckck ck ck ck ck ck ck ck ck ck ck ck k ck kk ethod Schwabisch Data output file Outdata hei schwabisch Data output format real4 First line w r t original master 1001 Last line w r t original master 2105 AA Ges a OOS abso is e TT eil E OAL Last pixel w r t original master 700 Multilookfactor azimuth direction YO Multilookfactor_range_direction Z Ellipsoid name a b WGS84 6 37814e 06 6 35675e 06 ARM AA Fe Fe FEE Kee Ne oe e RASCA Ke eo eoe ee ete oe eO eoe ee Ie ie a te ie Pee ee ee RAR A he he te ie he xt x End_slant2h _NORMAL KKEKKKKKKKKKKKKKKKKKKKK KKK KK KK KK KKK KKK KKK KKK KKK KKK Sk Sk Sk Ck Ck Sk Ck Pk AA KAR In the output files the height is stored number of lines etc multilooked unwrapped interferogram etc 32 3 Implementation 32 3 1 Method ambiguity This method yields first the heights of the line pixel and the position P x y z In this manner the geocoding can be done in the same step by converting P x y z
2. ck c 0c ck Ck ck kk Ck ck ck Ck ck ck ck ck ck ck ck Ck ck kc ck cock ck ck ck ck ck ck ck ck ck ck ck ck ck AAA RARA ck ck ck ck ko kk ko ko ck Sk kv kv ok ok End master timing NORMAL Ck ck ck ck ck ck ck Ck ck ck ck ck ck ck ck KKK kk ck Ck ck ck KKK KKK KKK KKK KKK KKK RA kk ko kc ko ko k kock AAA ko In the logfile some additional information is written such as estimated offset and correlation for each window position Estimated offsets equal to NaN are ignored during analysis 7 3 Implementation The master timing error is computed in three steps First for each window position the zero meaned master amplitude window is shifted over the zero meaned simulated amplitude and the correlation is obtained see equation D 24 by computing all pointwise products and dividing by the norms of the particular windows 26 Second we select the most frequently occurring line and pixel offset pair by searching for the highest fre quency of window offsets above the average correlation Finally we convert the obtained integer line and pixel offset to the master timing error using the following formulation lines 7 1 PRE 7 1 master_azimuth_timing_error pixels 2 RSR ve master_range_timing_error where PRF is the pulse repetition frequency and RSR is the one way resampling rate in Hz The algorithms which are used to compute correlation are the same as the ones used in the step COARSE CORR see th
3. matrix Type gt matrix Type gt matrix Type gt matrix Type matrix Type gt matrix Type gt matrix Type matrix Type matrix Type gt matrix Type void void void void void void void void void void uint uint lO uint IN uint pO uint pN const return number of pixels resize uint 11 uint pi setcolumn uint pixel const matrix lt Type gt amp COLUMN setcolumn uint pixel Type scalar setdata Type w sets 2w setdata uint li uint p1 const matrix lt Type gt amp A setdata window winin const matrix lt Type gt amp A window winA setdata const matrix lt Type gt amp A window winA setrow uint line const matrix lt Type gt amp LINE setrow uint line Type scalar showdata const show all data size const return nsize in matrix matrix lt Type gt sqr const matrix lt Type gt amp A friend matrix lt Type gt sqrt matrix lt Type gt sum const matrix lt Type gt 8A const matrix lt Type gt amp A int32 dim friend void writefile friend void wshift matrix lt Type gt 8A int32 n matrix Type gt matrix 157 Annex F Adding a module In this annex a description is given how to add a module to the Doris software First get general idea of the structure of the software source2html It is preferred to stay in same format as us General 1 Read input cards and parameters for your module 2 Add you
4. Figure 29 9 Buffering of complex interferogram in blocks for phase filtering For a block pixlo pixhi e g 0 15 the output equals fo ran overlap 23 pixlo overlap pixhi overlap 3 12 The number of output equals size 2overlap pixhi pixlo 1 2overlap 10 98 Chapter 30 UNWRAP In this chapter the processing of step UNWRAP is described This step is currently not implemented within the Doris software To obtain the unwrapped interferogram you should use another software for example one of the routines of Ghiglia and Pritt 1998 which can be obtained by ftp at ftp wiley com public sci_tech_med phase_unwrapping These software should not be considered public domain you ought to buy the book The slant to height con version and geocoding can only be done with an unwrapped interferogram Recently snaphu of Curtis Cheng was put in the public domain It is recommended you install this software as standalone executable and continue with Doris for geocoding afterwards METHOD snaphu can be used from within Doris But experience has to be gained how this software best performs Sometimes the coherence as computed by doris seems to contian NaNs not a number snaphu does not expect this and exits when this happens In Matlab the created coherence file can be easily corrected with e g q freadbk 9192 6687 coh 2577 f10at32 idxx isnan q idx where idxx 1 q idx 0 0001 Wieder ke epo V Xem Sato imcia
5. i a sinc rect 7 65 28 0 lt lr 1 1 lt 2 lt 2 2 lt lx 3 3 lt x 21 3 21 4 21 5 21 6 21 7 21 8 Chapter 22 FILTRANGE In this chapter the processing of step FILTRANGE is described This optional step filters the spectra in range direction of master and slave to reduce noise in the interferogram The noise reduction results from filtering out non overlapping parts of the spectrum This spectral non overlap in range between master and slave is caused by a slightly different viewing angle of both sensors The longer the perpendicular baseline the smaller the overlapping part Eventually a baseline of about 1100 m results in no overlap at all the critical baseline for ERS Assuming no local terrain slope A reduction of typically 10 20 in the number of residues can be achieved Method porbits filters based on the orbits perpendicular baseline for a constant given terrain slope Perform this step after coarse coregistration since the approximate overlap is used to filter both images The output images are cropped to this overlap To filter on the save side i e not to filter out too much use a negative terrain slope of e g 10 degrees This step is not recommended except perhaps to improve the coregistration polynomial for long baseline pairs After the resampling the range filtering then could be repeated on the original data with the adaptive algorithm Method adap
6. x x x x real8 c void tr2pix real8 rangetime real8 tri compli16 cr4toci2 complr4 x coarseporbit void coarsecorrel void coarsecorrelfft real4 corrfft matrix lt uint gt distributepoints void getoffset void finecoreg real4 coherencefft real4 coherencespace void coregpm matrix lt real4 gt getofffile matrix lt real4 gt cc4 matrix lt real4 gt cc6 matrix lt real4 gt ts6 matrix lt real4 gt ts8 matrix lt real4 gt ts16 matrix lt real4 gt rect matrix lt real4 gt tri void resample void rangefilter void rfilterblock void phasefilter matrix lt complr4 gt goldstein matrix lt real4 gt smooth matrix lt real4 gt smooth void spatialphasefilt matrix lt complr4 gt convbuffer void phasefilterspectral matrix lt complr4 gt spectralfilt void azimuthfilter matrix lt complr4 gt blockazifilt void slant2hschwabisch void slant2hambiguity void slant2hrodriguez void geocode void printcpu void inittest bool doinitwrite void initwrite void updatefile void getanswer bool readres void updateprocesscontrol void checkprocessing void checkrequest void fillcheckprocess void fillprocessed int32 filelines inline ONOMOSONSOEONEONONOSOMONSOSONRONMCOSOSONC 135 rangesamplingratex2 84 existed ioroutines c bool existed 85 removedatleader ioroutines c void removedatleader 86 filesize ioroutines c
7. y loat S12 ee If you use a standalone application to unwrap the interferogram you might have to mimic the output as described below so Doris can obtain the current filename and dimensions for the unwrapped interferogram from the interferogram result file 30 1 Input Cards For the snaphu program please also refer to their website and read the man page Doris is a wrapper for a system call to the executable snaphu Therefore a program called snaphu should be executable and in your path An input file for snaphu is created in the current directory You can rerun snaphu with a changed inputfile from the prompt if required but keep same output file name format as it is in the result file for the interferogram We assume to unwrap the complex interferogram mph format always for snaphu UW METHOD SNAPHU TREEFRAMON Select method for unwrapping For general users if they have installed the Snaphu program this will make a system call Other methods are not available fur public domain UW OUT FILE uint hgt Output filename of unwrapped interferogram 99 UW_OUT_FORMAT HGT REAL4 Output format of unwrapped interferogram UW SNAPHU MODE TOPO DEFO SMOOTH NOSTATCOSTS Output format of unwrapped interferogram Snaphu options t d s respectively Refer to snaphu manual for more information UW_SNAPHU_COH filename use specified file for correlation values This file must be registered and have the same number of lo
8. 1 add definition of new struct 2 readinput augment with new struct 3 add in big switch what to do if new step is requested General document what you did new keywords and arguments new process control flag e how does the result file has to end END filtphase NORMAL what strings in the result file are used later in the program e email to owner doris_users tudelft nl 160
9. B baseline R remark A author I h For more help type run h C 2 3 cpxfiddle With this utility one can fiddle about with binary complex cpx files of all kinds of formats though only pixel interleaved Make a cutout multilook scale exponent subsample mirror etc Now it also supports the 139 generation of SUNraster files for visualization of the phase of complex files smaller temp files required It is written in C using a template function Input is a complex file for example the output of Doris or SLC data It should be pixel interleaved i e RE IM RE IM RE IM Almost all binary pixel interleaved formats are supported Output is written to stdout channel normally your screen in ascii or binary float format The binary output should be redirected to a file or piped to a GMT program Ascii output can best be used to view a small cutout Output option 0 are normal the file as is the magnitude the phase the real or the imaginary part Further options are making a cutout multilooking subsampling and or mirroring in the vertical and horizontal plane Cpxfiddle does not handle band interleaved complex data column major order files nor non complex files However cpxfiddle might be tricked See also cpxfiddle h and the utility cox2ps C 2 4 SYNOPSIS l cpxfiddle w width f informat q output o outformat e exp s scale line L line p pixel P pixel
10. If such files are not removed by Doris for example after an error in the processing they can be savely removed by a rm command In the logfile additional information is written that is not sent to standard out or the result files such as statistical information on a least squares estimate For example we always inspect the correlation value for step COARSE CORR and the error between model and observations for reference phase computations in the log file PRINCIPLE OF RADAR COORDINATES and MULTILOOKING IN DORIS linelo linehi 1 2 3 4 5 6 7 8 9 10 11 12 13 14 eo e 09 09 9 9 master grid ee M L cu XM A complex interferogram 0 14 multilook 4 numelements on disk 3 complex interferogram refpha multilook 8 numelements on disk 1 Figure 1 3 Principle of multilooking 1 2 4 Inputtfile There is one ascii input file containing cards and zero or more parameters that controls the processor A card is the first word on a line in the input file The card and parameter s are delimited by blanks or tabs There are mandatory cards such as STOP at the end and optional cards because there are defaults for example for a filename for the output or cards like COMMENT In this document the mandatory CARDS are in sans serif bold face the optional CARDS are in normal style Default parameters are underlined and optional parameters are betwee
11. interfero coherence comp refphase subtr refphase comp refdem subtr refdem filtphase unwrap stant 2h geocoding dinsar NOTAUSEDAR End_process_control DESS SE SS SY SS O e SYS RSS SS SKIP SKIP The latter flag is reserved for future use As already mentioned after the process control flags the results of a successfully ran processing step are appended A section of the tail always starts with some lines like ck ckckckckckckck KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK ck ck ckck kk ckckck z _SIECUCE _COcigsS exa ole sis EXHXERAHARXAXE RL RA ko kc RARE RRE RA ERA ARA AA AE RAE kck ARA ARA ck k kk k kk kk After which the results for this processing step follow A section always ends with a statement shown below This End_ step NORMAL statement is important because the status of the process flag in the header is updated with it ok ck ck ck ck ck ck ck KK KK KA KK KARA KA ck ck ck AAA ck ck ck AAA ck ck ck kk AAA AAA AA ck ko ko ck ko Sk kv kx ko ko ko x End coarse orbits NORMAL Ck ck ckckckckckckckck ck ck ckckckckckckckckckckckckckckckckckckckckckckckck ck ck ck ck ck ck kck k ck ckckckck ck ck ck ck ck ck ck kc ck k ck kk kk Note that not all steps that are in the process control flags actually have to be implemented in Doris for the moment Unwrapping extra flags Since only one result section is allowed for every processing step it is not possible to r
12. known h to phi lambda If there is a trend in the height this has to be removed first e g by using tiepoints This means the computed phi and lambda matrices are not correct anymore With the baseline defined as in Annex D the following equations hold Bi T1 T2 32 1 By Bsin 6 a 32 2 B Bcos 0 a 32 3 Note the sign must be minus 4m i yr obj 32 4 Content of unwrapped interferogram never mind the phiR 4m 1 2 oR z PI OR 32 5 do 4rdBy 4a Bw XP Ee Geometric equation h H r1cos0 32 7 110 dh 97 sind 32 8 Height ambiguity dh dh d Arysn r sin 32 9 do d d Ar Bj 4nB cosO B sin0 Final equation to convert phase to height A ri sin 0 LLL 2 1 4v By cos 0 By und ipeo The procedure to compute the height is as follows note that computation is skipped if unwrapping went wrong indicated by NaN not a number in the unwrapped interferogram 1 for all lines a i1 compute B and B by computing h for middle pixel then compute baseline b i2 hcurrent 0 c for all pixels i j1 if value unwrapped phase equals NaN 999 then goto next pixel ii j2 compute 0 corresponding to h ii j3 hlast hcurrent compute hcurrent with formula v if hlast hcurrent gt k then goto j2 32 3 2 Method rodriguez Use same definitions as exact method see also Annex D See also Rodriguez and Martin 1992 This me
13. 82 Chapter 26 COMPREFDEM In this chapter the processing of step COMPREFDEM is described A DEM is radarcoded at the grid of the interferogram In step subtrrefdem chapter 27 it can be subtracted from the complex interferogram This step requires a DEM It can best be performed after the interferogram generation SRTM can be used to obtain the topography in 3 arcseconds Near global 90 m at the equator or 1 arcsecond USA only resolution The input DEM file must have the byte order of your platform in order to extract correct elevation value see CRD_IN FORMAT option for details The DEM should be in the WGS84 system same as the orbit ephemerides The Doris distribution contains the utility construct dem sh to download and prepare SRTM data see Section C 2 12 An alternative is the USGS gtopo30 DEM This is a global DEM with relatively low precision and gridspacing 30 seconds approximately 1 km at equator There are 33 tiles covering the globe in total requiring about 3GB of disk space For more information see for example http edc usgs gov products elevation gtopo30 gtopo30 html or Google for more in formation The input DEM has to be in an equiangular grid The format is short signed integers real4 floats or real8 doubles meters The DEM should be in WGS84 system same as orbit ephemerides The UpperL eft pixel of the DEM matrix on file should be the most North most West pixel i e the pixel with largest latitude
14. and combining this with 31 15 yields for the flat earth corrected phases 00B o An Bio 4m _ MD T EM 1 1 KSR ee a ue or Bio Ar g 14 feri eu 31 16 or for the phase ar caused by the deformation Ar B 1 dar P F9 31 17 LO This important equation shows how to obtain offset vectors from 3 SLC images i e by scaling the reference phase corrected unwrapped phase of the topo pair by the ratio of the perpendicular baselines to points on reference body and subtracting this from the phase of the defo pair This can thus be performed without the true values for 0 are required 31 3 1 Algorithm Input is the unwrapped topo interferogram corrected for flat earth Format is hat or real4 Not unwrapped thus indicated by NaN 999 real4 or if amplitude 0 hgt Defo interferogram is wrapped complex real4 mph specified in interferogram result file 1 Obtain orbit for topo slave result file Obtain filename dimension of unwrapped interferogram result file 2 Read in matrices appropriate size format etc per line Check if they exactly overlap 3 Compute B and B1 on a small grid 20x10 points over the image 4 Model the ratio of the perpendicular baselines by a 2D polynomial of degree 1 r l p a00 a10l a01p Give statistics for max error due to modelling It seems the ratio hardly changes over the image for ERS1 2 5 Compute wrapped deformation phase phase corrected
15. and other side of slave spectrum See also Fig 22 3 Note that the filter is mirrored matlab fliplr for master slave The SNR of the peak of a random spectrum sea probably is a little larger than 1 so threshold of 3 may not be large enough Do inverse FFT for filtered master slave which yields the filtered image in the space domain Take next part of master and slave e g 500 lines by next 128 range pixels until all lines are filtered In practice this is done in blocks These blocks are overlapping in lines because the averaging over lines means one cannot filter all lines in the block and not in range Parameters that can be adjusted are the FFTlength the moving average mean the SNR threshold The fftlength should be large enough to yield a good estimate of the local fringe frequency and small enough to contain a constant slope of the terrain The total number of fringes in range directorion can be easily estimated using the perpendicular baseline It is probably a good idea to add a card so an overlap in range between blocks can be used This avoids edge effects and increases the filtering of terrain near e g a lake since the SNR for peak detection will be higher for a number of blocks towards the noise This is not implemented yet See also Gatelli et al 1994 Geudtner 1996 Curlander and McDonough 1991 See also our matlab toolbox original spectrum master original spectrum slave 60 60 40 40 20 20 0
16. lt iostream gt has to be replaced by lt iostream h gt 3 The ios binary mode in file open is not required in IRIX 6 4 understand it s true for many unix variants simply removed it and then it works 4 In all declarations such as template matrix lt TYPE gt operator TYPE2 const matrix lt TYPE gt 8A const matrix lt TYPE gt amp B in matlib c the keyword template has to be replaced by keyword class This is also true for declaration for member functions such as template matrix lt TYPE gt correlate TYPE2 const matrix lt TYPE gt amp A const matrix lt TYPE gt 8B 5 During final linking for making doris executable an undefined symbol void matrix lt complex lt float gt gt conj was complained by the linker even though it was being generated by conditional comipling for creation of matrix library in matrixbaseclass overloaded This has mostly to do with order of linking the libraries avoided this by specializing the corresponding function void matrix lt TYPE gt conj as void matrix lt complex lt float gt gt conj in matrixbaseclass ovrloaded 6 get a lot of warning stating multiply defined malloc alloc oom malloc for different primitive types They can be ignored 7 A final error in ioroutine c while compiling in debug mode for line if compl4 1 1 complr4 1 1 0 stating more than one matches was reported Since it was only for debugging mode avoided this by commentin
17. slcimage cc slcimage cc r se Ge el Ge ee el readdata slcimage cc showdata slcimage h readvolume stepiroutines c readleader stepiroutines c readnull stepiroutines c readdat step routines c writeslc stepiroutines c unwraptreeframon unwrap c getorb utilities c convertgetorbout utilities c solve33 utilities c solve22 utilities c nextpow2 utilities c polyval utilities c polyval utilities c polyval utilities c polyvalid utilities c normalize utilities c normalize utilities c BBparBperp utilities c BBhBv utilities c Btemp utilities c BalphaBhBvBparBperpT10 utilities c iseven utilities h iseven utilities h iseven utilities h isodd utilities h isodd utilities h isodd utilities h ispower2 utilities h Ncoeffs utilities h degree utilities h remainder utilities h remainder utilities h sinc utilities h rect utilities h tri utilities h onedecimal utilities h onedecimal utilities h myrect utilities h myhamming utilities h interpbicubic utilities h normalize utilities h normalize utilities h void void void void checkcomprefdem checksubtrrefdem checkfiltrange checkdinsar void checkfiltphase void checkfiltazi inline void setunspecified char xs inline bool specified const char xs c void flatearth c void radarcodedem void slcimage fillslcimage const char file void slcimage updateslcimage matrix complr4 slcimage readdata inline void showdata const void rea
18. the filter for the slave is the inverse of that of the master The resulting spectra are shown in the bottom row The frequencies that did not overlap are filtered out yielding a better coherence between master and slave image The spectrum and filters depicted here are FFT shifted for clarity 41 e if the Doppler centroid frequencies do not vary per column use the same filter for all columns else compute the correct filter foreach column and use that First coarse coreg align then evaluate FDC polynomial e Take inverse 1DFFT in azimuth direction over the columns yielding the output The azimuth spectrum is also weighted for the antenna pattern fracsin fa foc foop fa fpc fpa 15 1 Where fDop 1505 Hz the Doppler bandwidth see Geudtner 1996 We did not de weight and re weight the spectrum for this This might be visible in figure 15 1 as a slightly asymmetric spectrum for master and slave We are not convinced that this re weighting can be performed without changing the signal However we believe that possible errors are small 42 Chapter 16 S FILTAZI In this chapter the processing of step S_FILTAZI is described Normally the step S FILTAZI is performed at the same time as M FILTAZI PROCESS M FILTAZI card in input file However we kept this two seperate steps to be able to only filter the slave images in a large stack all slaves coregistered on the same master image This has the advantag
19. treef_ramon ramon uw 100 100 ll kk ck kk kk S Sk Sk SS Sk Sk Sk S KK KK KK KK KK KK KK KK KK KARA RARA KKK KKKKKKKKKKKKKKKKKK End_unwrap _NORMAL KKEKKKKKKKKKKKKKKKKKKKKKKKKKKKK KKK KKK KKK KKK KKK KK KKK KKKKKKKKKKKKKKKKK If the Data output format is real4 then the output is assumed to be real4 unwrapped phase values If the unwrapping was not successful these pixels are set to 999 and ignored for slantrange to height conversion and differential insar 100 The Data output format of the unwrapped interferogram also can be hgt a band interleaved format ampli tude phase for SNAPHU For more details on the definition if unwrapping went wrong see Annex D 101 Chapter 31 DINSAR This chapter describes the processing step DINSAR which stands for 3 or 4 pass differential interferometry Three pass differential interferometry is described in Zebker et al 1994 It is a method to remove the to pographic induced phase from an interferogram containing topography deformation and atmosphere This module thus can also be used to study atmospheric effects in interferograms if no deformation is expected This step can be performed if an unwrapped topography interferogram topo pair and a complex deformation interferogram defo pair are present with a common master perform 2 seperate runs of Doris to achieve this The interferograms have to be corrected for the phase of the flatearth see ste
20. w r t original master i Last line w r t original master 29650 First pixel w r t original master 1 Last pixel w r t original master 4992 Multilcokfactor aczimubth direction 3 Multilookfactor range direction d Number of lines multilooked 29650 Number of pixels multilooked 4992 ck ck ck ck ck ck ke ck ck ck kk KKK ck ck KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK ko kk ko ko ko ko ko ko kv ko ko ko oko End sim amplitude NORMAL Ck ck ck ck ck Ck Ck ck ck kk KK KKK ck kk KKK KKK ck ck ck KKK KKK Ck ck ck kk kk ko kk kk ko kk kk ko kk Sk kc ko ko k ko ko ko kck ko The output in the logfile is more verbose specifying the results of the intermediate steps Also go over the standard out in case of problems with the SCREEN set to DEBUG level 23 Figure 6 1 The orginal multilooked amplitude of Figure 6 2 The simulated amplitude of the mas the master image ter image 6 3 Implementation The simulation of master amplitude is computed in three steps First the DEM is radarcoded to the coordinate systems of the master acquisition For each DEM point the master coordinates and the look angle 0 are computed and saved to temporary files Both the master coordinates and look angle are real valued Second the look angle and topographic heights are interpolated to the integer grid of master coordinates A linear interpolation based on a Delaunay triangulation is used The software package Triangle for the Delaun
21. 0 p 40 peak aetirmalitfl 150 0 50 100 150 5 0 filter6QSpectrunt Waster 300 filtered spectrum slave 60 100 40 50 20 0 0 0 50 100 150 0 50 100 150 Figure 22 2 Peak estimation in spectral domain of oversampled complex interferogram Non FFTshifted 70 22 3 3 Hamming filter The Hamming filter that optionally is used to de weight en re weight the spectrum of master and slave has the form W fr a as fls costar rect l 22 3 Where f is the frequency axis fs 2 df fs df df fs N f is the range sampling rate 18 96MHz and B is the bandwidth in range 15 55MHz a is a parameter controlling the amount of weighting fi llzli lt os wei 0 otherwise 224 Note rect not periodic filters 0 20 40 60 80 100 120 140 Figure 22 3 Spectral filtering windows inverse hamming boxcar rect and new hamming Note these are FFTshifted 71 detail of phase of interferogram before rangefiltering 20 40 60 80 1 00 120 140 detail of phase of interferogram after rangefiltering 10 15 20 25 j 30 35 40 20 40 60 80 100 120 140 Figure 22 4 Detail of interferogram with and without rangefiltering fftlength 128 nlmean 15 snrtresh old 5 The perpendicular baseline is about 200 m for this interferogram The fringes are clearly sharper after the filtering The number of residues for the interferogram was reduced by 2096 Subtraction of both interferograms yielded
22. 00 Initial offset for coarse coregistration The given offsets will shift master over simulated amplitude Use for debuging 0 0 by default Example input cards for this step C comment _ COMPUTE MASTER TIMING ERROR____ E c MTE_METHOD magfft computes faster than magspace MTE_METHOD magspace let audite MTE ACC TASAS only for magspace MTE_NWIN 2516 number of large windows MTE_WINSIZE 25 123 rectangular window MTE_INITOFE ONO initial offset 7 2 Output Description The process control flag at the start of the master result file is switched to 1 at successful exit m_timing dl Example of output of this step master result file Ck ck ck ck ck ck ck Ck ck KKK KKK KKK KKK KKK KKK KKK RARA ck kk kk ck ck ck kk ko kk kk ck kk ko kc ko ko k ko ko kc k ck ko Start master timing ck ck ck ck ck ck kc ck ck ck ck cock kk ck ck ck ck ck ck kc ck ck AAA RAR ck ck ok cock ck AAA AAA ck ck RAR AAA ck Ck ck sk ko Sk ko Mk kc ko ko ko Correlation method magspace 321 193 Number of correlation windows used JE o r 110 Estimated translation master w r t synthetic amplitude master dem Positive offsetL master image is to the bottom Positive offsetP master image is to the right Coarse correlation translation lines Coarse_correlation_translation_pixels 5 Master_azimuth_timing_error 0 000595265 sec Master_range_timing_error ESI eoo Ones
23. 001 Number of lines original 26183 Number of pixels original 4900 XN eoe oe RARAS AR Ie See EERE oe See eoe oe Seo eU ce eO eoe eae RAR AA SEEN e eS ee Ses Ie x End readfiles NORMAL KKKKKKKKKKKKKKKKKKKKKKKK KK KK KK KK KK KK KK ARA kckckckckckckckckckck ck ck ck ck ck ck ckckckckck k ck kk Note that Product specifies ASAR for ASAR which is used later 15 A number of these lines is only for your information The lines after Leader file are used for the further processing except Weighting identifiers These strings may not be altered We encountered a problem once when there was a blank in the UTC time instead of a zero but this should be fixed The logfile shows more details Also some information is echoed to the screen as INFO such as the Doppler centroid frequency evaluated at some ranges the corners of the images in latitude longitude etc Defaults for parameters slcimage cc Doris can still crash if not correct e g the approximate coordinates of the scene wavelength 0 0566660 Hi d eeu ES t rangel 5 5458 330 22 083 s one way default ERS2 Bye 1679 902 IV eal Clestapillic ies abw 1379 05 fi We cles JS ncs 18 9624680x2 0e6 ij zl certavit IRS 2 rbw 15552665 B clssusdb IS 3 3 Implementation Three basic data formats are supported CEOS ERS JERS ALOS RADARSAT N1 ENVISAT and COSAR XML TERRASAR X By specifying the M IN METHOD the correct reader is selec
24. 1 d 1 A0401 2 3 d 2 A290411 A02 4 5 6 d 3 A30A21 412403 7 8 9 10 77 Thus the number of coefficients unknowns equals s d 1 441 24 5 A polynomial of degree 5 normally is sufficient to model the reference phase for a full scene A lower degree might be selected for smaller images which also should increase the stability of the normalmatrix We would expect the higher order terms to be small because the polynomial describes a smooth long wave body ellipsoid To force the polynomial to be smooth one might consider always using a polynomial of degree 2 or 3 See also Schwabisch 1995 or Geudtner and Schwabisch 1996 24 1 Input Cards FE METHOD porbits Method selector for this step currently there is only one method FE DEGREE 5 degree of 2d polynomial FE NPOINTS 501 Number of points to compute reference phase for least squares estimation FE IN POS filename ascii file with positions master coord system in it where to compute the reference phase and then to model it by a polynomial Card FE NPOINTS is ignored if this card is specified This card can be used e g if it is desired to have the points on a grid including the edges Possibly one might even select points outside the grid though not smaller than 0 in order to avoid excessive fluctuations at the edge if a higher degree polynomial is used FE OUT FILE filename Card will be added in future for optional output of reference p
25. 180 180 Default is set to tile w020n90 DEM It is interpreted as max latitude min longitude in source SAM_IN_NODATA 9999 22 Identifier to ignore data in input DEM with this value Default is set to tile wO20n90 DEM SAM_OUT_FILE master sam Filename of the output simulated amplitude SAM_OUT_DEM filename Request optional debug output to float file of inout DEM per buffer cut to the inter ferogram window Info on these files is written as DEBUG Default is set to dem crop_sam raw Example input section E comment ___SIMAMP____ e SAM IN DE final WAna dem input DEM SAM IN FORMAT r4 real4 format SAM IN SIZE 3601 3601 extend of the DEM 1 p SAM IN DELTA 0 000833333 0 000833333 3 arcseconds SAM IN UL a 29 the center coordinates of the UL corner pixel SAM IN NODATA 32768 ignored value SAM_OUT_FILE Outdata 18226 sam synthetic amplitude SAM_OUT_DEM Outdata dem_sam raw cropped DEM for 6 2 Output Description At successful exit the process control flag is switched on sim_amplitude il The output in the master result file looks like Start sim amplitude Ck ck ck ck ck Ck Ck Sk ck kk KK KKK KKK KKK ck kk ck ck ck KKK KKK Ck ck ck kk kk ko kk kk ko kk kk ko kk Sk kc ko ko KKK KKK KK DEM source file final WAna dem Hm ee omes e D ENE ii AO PUED EME 252 Data_output_file Outdata 18226 sam Data_output_format real4 First line
26. 3 x 4 x Sun compilers Intel compilers and so on If you have problems installing Doris you can sent your questions to the email list of doris users To join this list follow the directions at our internet site http enterprise Ir tudelft nl doris Please don t forget to specify platform compiler versions etc B 1 Installation of Doris B 1 1 Installation of the Doris core After downloading the gzipped tarred archive of the Doris software v4 02the installation is best done with a Makefile assume you are familiar with make to compile code If you are not find someone who is 1 Create a directory for Doris e g mkdir opt doris v4 02 cd opt Doris v4 02 2 Download the archived Doris software via the download area of our webpages at http enterprise r tudelft nl doris doris v4 02 tar gz 3 Expand compressed files gzip d doris v4 02 tar gz This leaves a file doris v4 02 tar 4 Extract the files from the archive tar xvf doris v4 02 tar 124 Now sub directories are created bin src SARtools ENVISAT_TOOLS The doris source files are in src and bin The other two are utilities that need to be compiled separately Now we are ready to compile Doris cd to the src directory and read the README file for the latest information Create a Makefile by running the script configure in the sre directory csh configure or chmod 755 configure if it is not executable This should limit the editing in the
27. 4 1 Computation of coefficients o oo a D 4 2 Evaluation Of polynomials e vi D 5 SAR System parameters 149 D 5 1 AZIMU iso geo a a eee RW Oe we BO Ro e ec uid s 149 D 5 2 RAMOS a i ake oe AN 150 D 6 Doppler range and ellipsoid equations les 151 D 7 Orbitinterpolation moss 152 D 3 Format of the products cmo p ware Rr RO ee ev Ros Y Y x 153 E Matrix class 155 x ME En ee ELE MUNCIONS lt 2 45 aie ex a amp me Lr 155 F Adding a module 158 Edil on uri gt eye ob Bt os GE ew deem oh ch eee a EO a ee oe a N 158 E2 Adding a Step 26 0604 4s eb Re RRR ee ee RR Bee T e ew 159 Chapter 1 Introduction This user s manual guides you through the radar interferometric processing with the Doris software A user will be able to process ERS1 2 ENVISAT JERS RADARSAT ALOS and TERRASAR X data with the aid of this manual The data must be in Single Look Complex format SLC Doris is not a SAR processor i e you cannot process focus RAW radar data The manual also contains implementation specifications which can be helpful to a programmer for further development of Doris or to understand WARNING and ERROR messages 1 1 Overview of the InSAR Processing A high level description of the InSAR processing is shown in Figure 1 1 where a division in four blocks has been made Block depicts the preprocessing of the raw radar and orbit data to another format We won t be
28. First pixel w r t original master il Last pixel w r t original master 1000 ok oe oe cox eoe koe ok te ke ak tee Xe eh ose ao ie ee ak de Ko Xo cede OK ok ele KU X ect KC Ao Xo ko ko ak ake ko ko HE R End crop NORMAL Ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck Ck ck ck ck Ck ck ck Ck Ck Ck ck Ck Ck Ck Ck Ck CK CC CC AAA RRA RAR Sk Sk Sk I ke S ke ke kx kc k kc kc kc KKK If the SLC data is already on disk for example because the SAR processing was done this section will have to be simulated As well as the result from READFILES The format complex real4 is available Note that the byte order must be the same as the host platform order This means that if data is copied from big endian platforms they have to be swapped Use for example a dd command like dd if cdrom file slc of file slc conv swab or use gmt Zh for short Zf for float see man pages xyz2grd cdrom file slc Zh Sfile slc 21 Chapter 6 M_SIMAMP In this chapter the processing of step M_SIMAMP is described where the amlitude of the master image is simulated This step requires an external digital elevation model DEM such as SRTM Shuttle Radar Topography Mission elevation data This step can be performed optionally to simulate the master amplitude Eineder 2003 and should be done after the M_CROP step In the following step M_TIMING the synthetic amplitude can be used to compute the absolute timing error of the master acquisitio
29. RARA RARA I ke Sk ke ke kx kc ko kc kc kc KKK easy X m Y m Z m NUMBER OF DATAPOINTS 29 35360 000000 5161849 442 1645908 227 4678710 927 35361 000000 5166975 704 1645589 230 4673176 199 35362 000000 5172096 340 1645267 708 4667636 387 35363 000000 SS 1644943 664 4662091 497 35364 000000 BlLS2S20 FILS 1644617 098 4656541 536 35365 000000 5187424 437 1644288 011 4650986 509 35366 000000 MAS 1643956 405 4645426 423 35367 000000 SIRO NIGEL 1643622 281 4639861 283 35368 000000 S20 27 Ok 53 1643285 640 4634291 095 35369 000000 5207782 784 1642946 483 4628715 867 35370 000000 AMAS 204 1642604 811 MG 2S iL 3514 0S 35371 000000 21927 Sas 1642260 626 4617550 309 S537 2 5 OOOO 9222992 0102 MEAT 3 S28 AGES IZ 35373 000000 5228050 368 1641564 720 4606364 658 35374 000000 ASS LOS DS dL 1E 2719 AO 0 4600764 313 35375 000000 552 S8 3L 0 0 0 7 1640858 775 4595158 964 35376 000000 5243 TME 1640502 040 4589548 615 ISSN 7000000 5248226 814 1640142 800 45839337 213 35378 000000 5253256 642 16397 M054 4578312 945 18 35379 000000 5258280 744 1639416 805 4572687 636 35380 000000 52092 90RFAG5 1639050 053 4156 7051 58 35381 000000 ADS SINO 1638680 800 4561422 100 35382 000000 ASS a 1638309 047 4555781 886 ee oe Stoke she deo ERAN A RARA RA AREA AR AAN ARAS AR AR he te ie Pee RA o te Re TE Ne AN Ne Re AR AA RM Ke End_precise_orbits _NORMAL KEKKKKKKKKKKKKKKKKKKKK KK KKK KKK KKK KKK KKK KKK KKK KKK KK kk Ck Ck S
30. S x y 2M x y 2m mirror c file r rmin rmax B swap H bytes V b h elp ignorenan fullstat inputfile Dump content of complex binary file to stdout either as is magnitude phase real or imaginary part Input files can be almost any complex file though not yet band interleaved Output can be manipulated by multilooking subsampling mirroring scaling etc This program is useful for cropping and displaying in combination with e g GMT ImageMagick or xv Output format to stdout can be binary Careful only pipe or redirect this use cpxfiddle h amp more in csh or cpxfiddle h 2 gt 81 less in bash for more help C 2 4 cpx2ps With this utility postscript files from complex data and binary float data can be generated such as the output of Doris and SLC files Various input formats are supported Options are multilooking mirroring plotting the phase magnitude real and imag part etc It has become pretty big It calls requires cpxfiddle see C 2 3 and GMT subprograms grd2cpt grdimage psscale cpx2ps v2 1 FMH software Bert Kampes c 1999 2000 PROGRAM cpx2ps produce various encapsulated postscript code from complex data files SYNOPSIS cpx2ps w width f format cr4 q out mag e exp 1 0 s sc 1 0 1 TL alllines p 1 P width Mt 1 F1 1 T title c cptname gray z size 16 o epsfile G grdfil
31. Tro ble shooting 6 4 0 4 opus mem oom ee EO X he we eh ee de do d B 5 1 General problems eee B 5 2 Matrix class troubles eee ee B 5 3 Some notes on installation on SGI ll B 5 4 Some notes on installation on Linux X86 o e B 5 5 Some notes on installation on Window running Cygwin List offilesin archive 2222s List of routines description eA C Utilities C 1 C 2 C 3 Packages ura ae hb Owe ye tee whee Owe bea e ER ID ea A ee ee doe C 2 1 Installation of SARtools lll C 2 2 rur schlpt oe be See we so VOX XX RR eun a we ee Bod odes C23 eo C CC 2 4 CDX2DS 3o e Ruhe hah obese do D qose he ee croceo ELE E s edes 62 5 phasefilt doris lc m pores ROREM UV a ee Ree eg ogg C206 1lap ackz oye Rue PPM elei xe v ee Rb sue G27 SCPXMUNE PEE PC PPM C28 CPXdIV caba ARA mee ed ke A de ee eee RAG 6 29 CPC i a x x gue ete ete ea RO RR Re alin BBs WE a ea ple ien d Be C210 floatmull oxida soe A A RUE ERU RO Roe RAD AUR e Rt m Rd A Wrap s amp ek eem we RE RO ee ee wes uw dod des 2 12 constr ct dem str ux we BO eae e Roe ise Rm ee ee ee RR C 2 13 dOriS proCess Tesetsh a Completes for tesh users sls D Definitions D 1 D 2 D 3 D 4 GonstantS PC e BM Lee AA A E ULT Interferogram eie ai asb di ioi Si oe ae e a e d aaa ee P lynomialS sico A S eo e ade se BB U A e OR oe a D
32. Type operator const matrix Type amp A matrix lt Type gt 8 matrix Type operator Type scalar matrix lt Type gt operator const matrix lt Type gt amp A Type scalar matrix lt Type gt operator const matrix lt Type gt amp A const matrix Type amp B matrix lt Type gt 8 matrix Type operator Type scalar matrix lt Type gt amp matrix lt Type gt operator const matrix lt Type gt 8A friend ostream amp operator lt lt ostream amp file const matrix lt Type gt amp A matrix lt Type gt amp matrix lt Type gt operator const matrix lt Type gt amp A bool matrix lt Type gt operator Type scalar const bool matrix Type operator const matrix lt Type gt amp A const friend istream amp operator gt gt istream amp file matrix lt Type gt amp A Type matrix lt Type gt operator uint line const Type amp matrix lt Type gt operator uint line uint pixel const matrix lt Type gt matrix lt Type gt operator window win const 156 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 operator pixels readfile resize setcolumn setcolumn setdata setdata setdata setdata setrow setrow showdata size sqr sqrt sum writefile wshift matrix matrix Type matrix Type operator uint matrix lt Type gt pixels const friend void readfile matrix Type amp Result const char file matrix Type
33. a random phase so no structural effect of range filter implemenation is suspected 72 Chapter 23 INTERFERO In this chapter the processing of step INTERFERO is described In this step the following is computed The complex interferogram is computed with or without subtraction of the reference phase The reference phase is subtracted if there is a 2d polynomial in the products result file result of step FLATEARTH It is not subtracted if this is not in the result file or if the number of coefficients is set to 0 The complex interferogram minus reference phase is defined as I 2M S R 23 1 Where denotes the complex conjugated denotes a pointwise multiplication lis the complex interferogram M is the complex master image S is the complex resampled slave image Ris the complex amplitude 1 reference phase The phase image of the complex interferogram minus reference phase is defined as arct no Juaag treat 23 2 Where arctan is the four quadrant arc tangent o is the phase image lis the complex interferogram This is identical to 0 w 6s Or 23 3 Multilooking can be performed to reduce noise Usually a ratio of line pixel 2 5 1 between the factors is chosen to obtain more or less square pixels 20x20m for factors 5 and 1 The resolution decreases of course if multilooking is applied 73 1 1 i 1 hh i We il 1 H In Il
34. and this step yielding new coefficients The polynomial coefficient can then be added to form a new model with which the slave can again be resampled This has not been tested 20 1 Input Cards CPM_THRESHOLD 0 4 Threshold for correlation value to use estimated offset of step FINE in estimation of polynomial coefficients This depends on the size of the window during FINE Estimated coherence using small windows are more biased towards 1 0 so a higher threshold is better For window size 64 64 a threshold 0 2 seems OK The plotoffsets script can be used from the prompt to figure out a good threshold value CPM_DEGREE 1 Degree of 2d polynomial See annex for definition of degree Degree 2 is adviced CPM DUMP OFF ON Dump computed model to files in float format Filename for azimuth model is offse tazi_ l p r4 where number of lines pixels are substituted Filename for range is sim ilar Content of file is evaluated model in master system via INFO the dimensions are also echoed CPM PLOT NOBG BG 57 Call gmt script plotcpm to plot results and to view with gv An example of the plots is given above The argument NOBG prevents a call to cpxfiddle to generate a magnitude background while BG does call cpxfiddle See the script plotcom and the c program cpxfiddle for more information cpxfiddle can be downloaded from Doris internet pages The command is echoed to stdout as INFO which can be repeated outside Doris CPM_WEIGHT
35. and unzipped as from version 2 5 This step introduces a section in the result file where the ephemerides are placed and it deletes the ephemerides from the SLC leader file obtained by the processing step M READFILES if there was such a section The ephemerides x y z span the time 4 seconds before the first line and 4 seconds after the last line by default The time is and should be in seconds of day The time interval is 1 second by default Natural cubic splines are used for interpolation and the boundary conditions may affect the interpolation if only a few datapoints are used e g 5 points with a time interval of 30 seconds We nowadays use a time interval of 30 second and approximtely 21 points This implies that only spline degree 3 piecewice polynomial is used for the whole image which gives better results for e g reference phase computation The interpolation errors in Doris are probably due to interpolation of interpolated values of getorb which output format is in 3 digits If you want to use other ephemerides you can simply insert them in the result file in the format described in section 4 2 You will have to correct the number of POINTS in the result file Note that the orbit system is WGS84 only 4 4 Input Cards M ORBDIR directory name the tar archive directory name for the Delft Orbital Data Records M ORB INTERVAL 1 Time in seconds between ephemerides M ORB EXTRATIME 3 17 Time in seconds before first a
36. be viewed with any standard editor For this the run script utility can also be used customized by setting the PAGER and EDITOR environment variables The data output are binary files which can be visualized with standard software on your system e g Khoros Matlab tip try spinmap etc The utility cpxfiddle c in SARtools archive can amongst other things generate SUNraster files of the phase magnitude and magnitude phase with a colormap of your liking use e g xv or display to view print This seems to be one of the fastest ways to see the results See also the PREVIEW card The utility cpx2ps csh script can generate postscript code from complex files for the magnitude or phase use gv or ghostview to view print This utility uses the GMT package Tip generate 2 m postscript files with Z option and view this enlarged interferogram with gv If someone want to develop a Motif Lesstif X windows application that would be very welcome We have experimented with these things but no time to implement it in a robust way A tip might be to use X utilities like xmag or xlupe to zoom in on the visualized files 128 B 5 Trouble shooting B 5 1 General problems In Tables B 1 and B 2 a number of possible problems and their solutions is given Table B 1 Troubleshooting 1 Problem Solution Function FFT not found while linking You have no Veclib library Change define statements in Makefile Run the install sc
37. concerned with the raw orbit data but it is included in the flow chart for completeness The Delft precise orbits are used for ERS1 2 and Envisat obtained via the getorb package see e g Scharroo and Visser 1998 The second block consists of the co registration where the slave image is aligned with the master image and of the computation of the reference phase of the ellipsoid In block Ill the interferometric products complex phase image and coherence map are computed Finally in block IV the endproducts e g a DEM or a deformation map are computed Figure 1 1 Coarse flow of the interferometric processing of SAR images The processing steps that are implemented in the Doris software are listed in the table below Note that UNWRAP is NOT directly implemented in Doris To run each step these names must be used as arguments for the PROCESS card as explained in Chapter 2 A specific algorithm module method can be selected with the cards that are special to this step see the Chapters 3 to 33 PROCESS Chapter Description Description amp amp O M READFILES 3 Read the processing parameters from the SLC files for the eS maser mage ee M_PORBITS 4 Retrieve the precise Delft orbital data records with the getorb ee eager M_CROP 5 Write the SLC data from paf format to disk in the raw pixel interleaved 2b 2b complex short integer format M_SIMAMP 6 M_TIMING 7 Estimation of ti
38. const matrix lt Type gt amp A Type min const matrix lt Type gt amp A Type min const matrix lt Type gt 8A uint amp line uint amp pixel matrix lt Type gt multilook const matrix lt Type gt amp A uint factorL uint factorP void matrix Type gt mypow Type s friend void myswap matrix lt Type gt amp A matrix lt Type gt 8B bool matrix lt Type gt operator Type scalar const bool matrix lt Type gt operator const matrix lt Type gt amp A const matrix lt Type gt operator const matrix lt Type gt amp A const matrix lt Type gt 8 B matrix lt Type gt operator const matrix lt Type gt 8 A Type scalar matrix lt Type gt operator Type scalar const matrix lt Type gt 8A matrix lt Type gt 8 matrix lt Type gt operator x Type scalar matrix lt Type gt amp matrix lt Type gt operator const matrix lt Type gt 8A matrix lt Type gt operator const matrix lt Type gt amp A const matrix lt Type gt 8 B matrix lt Type gt operator const matrix lt Type gt amp A Type scalar matrix lt Type gt amp matrix lt Type gt operator const matrix lt Type gt amp A matrix lt Type gt amp matrix lt Type gt operator Type scalar matrix lt Type gt operator const matrix lt Type gt amp A const matrix Type amp B matrix lt Type gt operator const matrix lt Type gt amp A Type scalar matrix lt Type gt operator const matrix lt Type gt 8 A matrix lt Type gt amp matrix
39. for topography with formula 31 17 using the modeled ratio r Actually compute it complex c cos r 9 i sin rij 6 Set not unwrapped regions to 0 0 7 Write output file complex real4 mph format if a problem with unwrapping occurred write 0 0 If requested also write the scaled unwrapped interferogram 106 Figure 31 4 Geometric configuration for 3 pass differential insar The orbits go into the paper All angles are defined counterclockwise The terrain element P corresponding to the radar coordinate l p is located at a height h above the ellipsoid The perpendicular baseline required for this method is the one for points P located on the reference ellipsoid h 0 00 the change in 0 since P is on a height due to a 5 km height difference is approximately 1 107 Chapter 32 SLANT2H In this chapter the processing of step SLANT2H is described In this step in principle the heights in the radar coded system are computed However with the exact method the geocoding can be done in the same step The results of the three implemented methods are different so a comparison has been made Processing is in buffers for all methods while it is possible just to do it line by line In case a polynomial has to be evaluated rodriguez method it is more efficient to have a buffer 32 1 Input Cards S2H METHOD ambiguity schwabisch rodriguez Method selector ambiguity geocodes as well u
40. however the step COARSECORR has to be run in order to get the coarse offsets within a few pixels FINE requires the initial estimated offset within a few pixels The offset is defined in such a way that for a point P in the master with coordinates Pm line pixel and the same point in the slave image with slave system coordinates Ps line pixel it holds Ps l p Pm l p offset J p 13 1 13 1 Input Cards There are no input cards for this step i e the parameters orbits are read from the master and slave result file 13 2 Output Description Since this normally is the first step that is not specific to master nor slave a products result file is created The process control flag at the start of this file is switched to 1 at successful exit coarse_orbits I Example of output of this step KEK ck coke ck ke KKK KEK RA EA RA KKK KKK ARA AAA AAA EA kok kockck ck ck kockok ck kk kok k kk x Start coarse coregistration based on orbits kckCckCckckckckckckckckckckck ck ckckckckckckckckckckckckckckck ck ckckck ck ckckckckckckckockckckckckckckckckck ck ck ck ck ck ck ckckck ck ck k ck kk Some info for pixel 3037 590 not used Bperp m 36 2 Bpar mii DDD Bh m 41 8 Bv m SEA B m 43 3 alpha deg 15 4 baseline orientation theta deg 18 look angle Height_amb m 20208 Btemp days i 34 Estimated translation slave w r t master Coarse_orbits_translation_lines BAG Coarse_orbits_translati
41. i ff LOGFILE log out M_RESFILE master out d S RESFILE slave out d I RESFILE interferogram out E AN SKIP SKIP more cards specific to step specified by see next chapters for details on these cards SKIP SKIP STOP Note that is not a delimiter for comments text after the last expected parameter is simply ignored 13 Chapter 3 M_READFILES In this chapter the processing of step M READFILES is described It can be selected by a PROCESS M READFILES line in the input file This is the first step if the ERS1 2 SLC images are processed The SLC leader volume and header of the data file are read and relevant parameters are written to the master result file specified by the general card M RESFILE These parameters are used in the further processing Currently ERS1 2 SLC and ENVISAT SLC files can be read If the output of this step is mimicked Doris can be tricked to process the other steps The sole purpose of this step is to create result file where relevant parameters are stored PRF wavelength etc also see the example in the next section 3 1 Input Cards M IN METHOD ERS ASAR ENVISAT RSAT RADARSAT ATLANTIS JERS ALOS TSX TERRASAR X Method selector to read ERS ENVISAT RADARSAT JERS ALOS or TERRASAR X header Note that both master and slave need to be acquired by the same sensor in principle JERS simply uses ERS programs ATLANTIS sar processor uses the ceos reader for RSAT and
42. ie aoa ar ai aO e e 21 3 Implementation aoaaa 0000002 21 3 1 Output formals ios a pre sanie O a a i i m urs 21 3 2 Interpolation Kernels 31 31 32 32 33 33 34 34 34 35 36 36 37 37 37 38 39 39 40 40 22 FILTRANGE 22 1 Input Cards 22 2 QUIDULIDESCHPHOM gt a oee m zee Rae Rom Bea Ee arse wee e een RUBRO De UR HR UR UB Re we ae 22 3 Implementation 22 3 1 Method porbits es 22 3 2 Method adaptive ry i e e ara ee ee ee Sa aes 22 3 3 Hamming filter 3c RR GG TR PR ERA pe ee he ee 23 INTERFERO 23 0 4 NEW 23 1 Input Cards 23 2 Output Description rs 23 3 Implementation 24 COMPREFPHA 24 1 Input Cards 24 2 Output Description x ocu Roe mum ec Y yp e ue mw e RO OY UE Y Rom x0 SR 25 SUBTRREFPHA 25 1 Input Cards 25 2 Output Description 2 2 2 08 o Roh eR ae UR EE Eod dre ee nds 26 COMPREFDEM 26 1 Input Cards 26 2 Output Description uos zoe iom a oO PUR Rok e ve edem roe FCR AR Rm Role S 26 3 Implementation 27 SUBTRREFDEM 27 1 Input Cards 27 2 Output Description uoo RUP E He irr m Re Uy Yep Sube qw a 28 COHERENCE 28 1 Input Cards 28 2 Output Description o 28 3 Implementation 29 FILTPHASE 29 1 Input Cards 29 2 Output Description e o uoo Be Ree ook mom eee cede A 29 3 Implementation 29 3 1 spatialeonv es 29 3 2 spectral 29 3 3 goldstein 30 UNWRAP 30 1 Inp
43. interfer ogram is related to the master this card has not a large effect STOP After this card the input file is no longer interpreted This card is mandatory 2 2 Example General Input Cards C KKKKK KKK KKK KKK KKK KKK ck ck ck ck ck ck ck ck ck ck ck ck ck ck ckck ck ck ck ck ck ck ck ck KKK KKK KKK kck ck ck ck ck ck ck AAA 12 2000 Monday x level of output to standard out MB overwrite existing files non interactive prevents copy of this file to log read parameters obtain precise orbits crop data to internal format log file parameter file parameter file parameter file ONLY PROCESS cards c Doris inputfile generated by run at Nov 27 C ckck ck ck ck Ck ck Ck ck ck ck ck ck ck Ck ock ck ck ck ck ock KAR ck ck ck ckock koc ock ko ck ck ckock ck ck ck kk ok ck ck ck ko ck ck ko ck ko ko Sk ok Mk Sk Mk kx Sk ko kx ok Cy c x Filename Inputfiles input s initial Ch Puro uc Doris User ic Master 23195 s Slave 05512 c Baseline JL sim c Remarks Test s2h routine exact Q C ckck ck ck ck kk ck Ck ck ck 0k ck Ck ck Ck ock ck ck ck ck ock oko ckock ck ock AA KA KARA KARA ck oko ck ck ko ck ko kx Sk ck Mk sk Mk kx Sk ko ok E comment _ general options___ e SCREEN debug Hi MEMORY 150 eh OVERWRITE BATCH ii c LISTINPUT OFF E PROCESS m_readfiles fd PROCESS m porbits Mik PROCESS m_crop Hil E Hf E comment _ the general io files Hf if E
44. is done with 2 spaces with the curly braces as shown below if expression action action2 You can display information depending on the value of the variable displevel with the functions DEBUG char ONE27 INFO char ONE27 PROGRESS char ONE27 WARNING char ONE27 ERROR char ONE27 F2 Adding a Step Adding a new step is not intended to be necessary The only thing that needs to be added are modules methods in pre defined steps However we will explain what you will have to change if you want to add a new step In file readinput h 1 you will have to add a const for the new step which is later stored in the process array 2 also a struct has to be made to store the variables of this new step method selector output file name window sizes etc 3 the prototype of the function readinput should be augmented with this new struct In the file readinput c 1 function readinput augment with new struct 2 only process card define new keyword for this step 3 add reading of parameters into defined inputstruct by new keywords In file ioroutines c only minor adding 159 1 routine doinitwrite add new step routine initwrite process control routine updateprocesscontrol check for string routine checkprocesscontrol check for string routine fillcheckprocesscontrol check for string D a RF wW routine fillprocessed check for string In file processor c
45. master 2105 iube Pi o lo ci a sie 50m Last pixel w r t original master 700 Ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck Ck Ck ck Ck Ck Ck Ck Ck C C C C CC C CC CC CK CIC RARA RRA KARA ke Sk I Sk S e AAA k kc kc kc kckok End resample NORMAL Ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck Ck ck ck ck Ck ck ck ck Ck Ck Ck Ck Ck Ck Ck C C CC CK CIC CK CC CK CK SS Sk Sk Sk ke Sk S ke S S kk kx ko k kc kc KKK 62 Soul overlap pixels Figure 21 1 Definition of the overlap between master and slave Note that the line and pixel numbers are given in the master coordinate system because the slave is interpo lated to that grid now 21 3 Implementation The overlap between slave and master is computed as indicated in figure 21 1 Interpolation is done with a kernel function such as e g a truncated sinc function First a look up table is computed for the selected interpolation kernel This table evaluates the kernel every 1 INTERVAL 0 05 positions which should be accurate enough Interpolation is independent for azimuth and range direction For all points in the overlap between master and slave the co registration polynomial is evaluated If the total image does not fit in the memory processing is done in buffers The azimuth spectrum of the complex SLC data is not centered around zero in general The location of the peak in the spectrum is given by a polyn
46. of interpolation kernels for SAR interferometry IEEE Transactions on Geoscience and Remote Sensing 37 1 318 321 Nitti et al 2008 Nitti D O Hanssen R F Refice A Bovenga F Milillo G and Nutricato R 2008 Evaluation of DEM assisted SAR coregistration In SP E Europe Remote Sensing Proceedings 15 18 September 2008 Cardiff United Kingdom pages 1 14 Rodriguez and Martin 1992 Rodriguez E and Martin J M 1992 Theory and design of interferometric synthetic aperture radars EE Proceedings F 139 2 147 159 Samson 1996 Samson J 1996 Coregistration in SAR interferometry Master s thesis Faculty of Geodetic Engineering Delft University of Technology Scharroo and Visser 1998 Scharroo R and Visser P 1998 Precise orbit determination and gravity field improvement for the ERS satellites Journal of Geophysical Research 103 C4 8113 8127 Schwabisch 1995 Schwabisch M 1995 Die SAR Interferometrie zur Erzeugung digitaler Gelandemod elle Forschungsbericht 95 25 Deutsche Forschungsanstalt f r Luft un Raumfahrt Oberpfaffenhofen 120 Shewchuk 1996 Shewchuk J R 1996 Triangle Engineering a 2D Quality Mesh Generator and Delaunay Triangulator In Lin M C and Manocha D editors Applied Computational Geometry Towards Geometric Engineering volume 1148 of Lecture Notes in Computer Science pages 203 222 Springer Verlag From the First ACM Workshop on Applied Computational Geo
47. order to trick Doris to use other parameters then the ones that result from a previous processing step simply edit the result file e g in order to coregister complex interferograms For the master output file the header with the information and the process control flags looks like MASTER RESULTFILE Master Tes Created by InSAR Processor Doris Delft o o Radar Interferometric Software Version Version 4 01 19 DEC 2008 optimal FFTW library used VECLIB library not used LAPACK library not used Compiled at Des 19 2003 178260852 By GNU gcc 4 1 4 File creation at eat Dee 19 SAO SZ 0 Delft Institute of Earth Observation and Space Systems Delft University of Technology http enterprise lr tudelft nl doris GE J9rcoxeess XE readfiles precise orbits uo pr sim_amplitude master_timing oversample calle alzado filt_range NOT USED End process control ey xe esr Sr xem KS Yes TS ei The last flag NOT_USED is reserved for future use In the slave result file the following line and processing step is extra resample 0 and following lines are missing sim amplitude 0 master_timing 0 The products result file is build up in the same manner The process control flags in the header of the products result file are SEA Ea Sig autem Oc SS cB coarse orbits coarse correl fine coreg timing error dem assist comp coregpm
48. remain some residual effects A next version of Doris will include an option to first correlate the radarcoded DEM with the interferogram to find an additional offset The radarcoded reference DEM is then shifted on multilooked pixel level before subtraction 88 Chapter 28 COHERENCE In this chapter the processing of step COHERENCE is described In this step the following is computed The complex coherence image is computed with or without subtrac tion of the reference phase and the radarcoded DEM phase The reference phase is subtracted if there is a 2d polynomial in the products result file result of step FLATEARTH It is not subtracted if this is not in the result file or if the number of coefficients is set to 0 For a proper estimation of the coherence map especially over mountainous areas the radarcoded DEM phase should also be subtracted In this case the result of step COMPREFDEM is required in the products result file This is the default procedure since v4 01 The complex coherence image between two images is defined as 7 E M S TU JAM M ElS 5 Where Ef is the expectation is the complex conjugated Ye is the complex coherence Mis the complex master image S is the complex slave image possibly minus complex reference phase and DEM phase S S R 28 1 The coherence is defined by y and its estimator as Nico MS x Xo MiMi 2 28 2 a Si Si Multilooking can be perf
49. resampling is required since the local fringe frequency is estimated from the interferogram This fringe frequency is directly related to the spectral shift in range direction Note this shift is not a shift but different frequencies are mapped on places with this shift The algorithm generally works as follows Take part of master and slave e g 500 lines by 128 range pixels Oversample master and slave and generate complex interferogram e Take FFT over range for all lines of complex interferogram e Take power If requested weight this powerspectrum with auto convolution of 2 rect functions with appropriate bandwidth Actually perhaps the spectrum should also be weighted with autoconvolution of Hamming but since am not sure that this has a big impact on real data this is not done e Take moving average over the lines of the power FFT s for noise suppression kind of periodogram This was better implemented as a convolution with a block function e g 9 x 128 Estimate peak per line in oversampled averaged powerspectrum of complex interferogram Estimate SNR fitlength peak rest 69 This peak is directly related to overlap of spectra for this part of this line See also Fig 22 2 Af Ftringe If SNR is above threshold input of user e g 3 remove appropriate parts of spectra of master slave Optionally compute inverse hamming window and new hamming window and rect window to filter one side of master spectrum
50. run Doris but it is highly recommended to use We choose to only use freely available packages C 1 Packages The program getorb for automatic retrieval of the Delft precise orbits for ERS1 2 See our web pages for a link The package gv recommended or ghostview to view postscript files generated by plotcpm and plotoffsets and cpx2ps The GMT generic mapping tools Highly recommended We generally generate postscripts for visualizing the phase and amplitude with this program See also cpxfiddle and cpx2ps C 2 Tools We have developed some utilities for running Doris and some display tools based on GMT C 2 1 Installation of SARtools After downloading the gzipped tarred archive of the SARtools the installation is best done with a Makefile assume you are familiar with make to compile code If you are not find someone who is See also annex B 1 Use a distinct directory to put the files e g mkdir opt doris_v4 02 cd opt Doris_v4 02 2 Download the archived SARtools via the download area of our webpages at http enterprise r tudelft nl doris SARtools tar gz 3 Expand compresed files gzip d SARtools tar gz This leaves a file SARtools tar 138 4 Extract the files from the archive tar xvf SARtools tar Now a subdirectory has been created SARtools Now the compilation of the utilities can start 1 Compile all utilities and install the executables cd SARtools inspect edit the Makefile set INS
51. step COARSECORR These computations are described in that chapter That step still has to be performed because the FINE step requires accurate initial estimates of offsets The AccL and AccP cards define the size of the searchwindow 2 AccL x 2 AccP around the initial offsets to interpolate a maximum The oversampling is done as follows 1 Transformation to spectral domain of searchwindow with correlation values at pixel level 2 Padd with zeros half the last term 3 Inverse transform 4 Find maximum in space domain this corresponds to estimated offsetvector Note that this way of computing is exactly the same if you first interpolate the signal and compute all correla tions and find the maximum or that you first compute at pixel level and interpolate the correlation values 47 Chapter 18 RELTIMING This chapter describes the processing step RELTIMING The relative timing error between the master and slave is computed using the precise orbits and the result of the fine coregistration Therefore this step can only be run after the steps COARSE_ORBIT and FINE The timing error in azimuth as well as in range direction is estimated This step is of influence for the coregistration using a DEM DEMASSIST step For this coregistration the master and slave images should be aligned as good as possible with the DEM The master is aligned with respect to the DEM using a simulated amplitude image M TIMING step This could also h
52. to stdout Consider using NOPLOT options in Doris input Read the help by typing plotoffsets at the prompt View the script in an editor try to detect the error adding x to the first line of the script should echo all commands before execution Run it without the background option What awk are you using We use gawk actually but it should be standard POSIX What is the correct syntax for tail n 6 or n 6 or 6 the tail commands have been removed from the new version This should not happen but sometimes these files are not re moved If Doris exits normally these files can be savely removed If Doris exits with an error sometimes these files can be used to repair the result files This ascii file is used to generate plots with GMT Try to find out with the debug version of Doris where it exactly goes wrong Some cards with optional ON OFF options need one of these specified on some systems 130 B 5 2 Matrix class troubles Some users reported having problems compiling a template class which was implemented in different files hope that by explaining my problems might help someone else For the aCC compiler the flag inst_implicitinclude had to be used This flag include the file matrixklasse cc note cc automatically Furthermore for making an archive library the member functions had to be instantiated explicitly see file matlib c for how did that The friend functions also had to be instantia
53. was performed and a histogram has been made of the correlation values The bin width equals 0 05 It can be clearly seen that the histogram of the filtered data is located to the right e is better 22 2 Output Description The process control flags in the result files for the master and the slave are switched on filt_range 1 In the filtrange section the filename and format is described after the filtering Master and slave are cut out so they exactly fit on each other KARA coke AR Hee Ke eoe Kee RER Tee A te Ke eoe ARA te Ke eoe GS ese RAR te Fee Ke ERREKA tee KO ae te de te ak ae tee Hee notarea lenrange KEKKKKKKKKKKKKKKKKKKKKKKKK KKK KKK KKK KKK KKK KKK KKK KKK kk Ck Ck Sk Ck Pk A A A A A ko ARA ethod adaptive Data output file utdata 23070 ee eres Data output format complex real4 mese laine sient Origa mal Emas EEr zT Last line w r t original master 6000 First pixel w r t original master 100 Last pixel w r t original master 3000 KKK ko kock ko ck c kokck c kockok ck E ARK ERE RARA ARE kokok kok kk kok ko ok ok ko k kok kk kokokck k kokokok kk ko ok k Kok End filt range NORMAL ck ck Ck ck ck ck ckckckckckckckck ckckckck ckckckckckckckckckckckckckckckckckckckck ckckckckckckckckckckckckck ck ck ck ck ck ck ck ckck ck ck k ck kk Figure 22 1 shows the improvement in the correlation for method porbits A histogram is made of the corre lation values of the FINE coregistration with and without range fi
54. 01161 Jul 17 47 2000 Source new ioroutines c r r 413 22 5289 Jul 17 47 2000 Source new ioroutines h r r r 413 22 29980 Jul r r r 413 22 70878 Jul r r r 413 22 57993 Jul 17 47 2000 Source_new matdoc txt 17 47 2000 Source_new matrixbk cc 17 47 2000 Source_new matrixspecs c r r 413 22 51646 Jul 17 47 2000 Source_new processor c r r 413 22 56881 Jul 17 47 2000 Source_new products c r r 413 22 2587 Jul 17 47 2000 Source new products h r r 413 22 130910 Jul 17 47 2000 Source new readinput c r r 413 22 15893 Jul 17 47 2000 Source new readinput h r r 413 22 42007 Jul 17 47 2000 Source_new referencephase c r r 413 22 3162 Jul 17 47 2000 Source_new referencephase h r r 413 22 1106 Jul 17 47 2000 Source new refsystems h r r 413 22 72665 Jul 17 47 2000 Source_new step1routines c r r 413 22 1121 Jul 17 47 2000 Source_new steptroutines h r r 413 22 11938 Jul 17 47 2000 Source new unwrap c r r 413 22 844 Jul 17 47 2000 Source new unwrap h r r 413 22 66820 Jul 17 47 2000 Source_new utilities c r r 413 22 9630 Jul 17 47 2000 Source_new utilities h r xr x 413 22 21633 Jul r xr xr x 413 22 13534 Jul r xr xr x 413 22 13093 Jul rwxr xr x 413 22 29621 Jul r xr x 413 22 548 Jul 17 47 2000 Bin helpdoris 17 47 2000 Bin plotcpm 17 47 2000 Bin plotoffsets 17 47 2000 Bin run 17 47 2000 Bin viewanddel AAA SAS
55. 2 bit and 64 bit platforms Processing of TerraSAR X data Improvement of the master slave overlap calculation The new algorithm should prevent segmentation faults as sometimes experienced in the past Improvement and speed up of COMPREFDEM see Chapter 26 based on a two step approach First the complete DEM is radarcoded and written to a file Second the radarcoded DEM is interpolated to the master image geometry using Delaunay triangulation The overlapping buffers as in previous versions of Doris are prevented Simulation of master amplitude image based on DEM M_SIMAMP see Chapter 6 Estimation of master timing error based on correlation between master amplitude and simulated master amplitude M_TIMING see Chapter 7 Estimation of relative timing error between master and slave image RELTIMING see Chapter 18 DEM assisted coregistration DEMASSIST see Chapter 19 For optimal performance apply after M_SIMAMP M_TIMING and RELTIMING Coherence estimation after DEM subtraction default see Chapter 28 Utility construct dem sh to download merge and fill voids of SRTM data see Section C 2 12 122 Utility doris process reset sh to reset and clean up the processing entries in Doris result files allowing to easily re run a step or multiple steps see Section C 2 13 Option to output height to phase conversion factors H2PH based on flat earth SUBTRREFPHA or DEM COMPREFDEM During COARSECORR and FINE coregistra
56. 29 Example input cards for this step 48 comment RELATIVE TIMING ERROR e RTE_THRESHOLD 0 4 RTE_MAXITER 10000 RTE_KALPHA 1 97 18 2 Output Description The process control flag at the start of the products result file is switched to 1 at successful exit timing_error 1 Example of output of this step products result file KKKKKKKKKKKKKKKKKKKKKKKKKK KKK KKK KKK KKK KKK kk Sk Ck Ck Ck Ck Ck Ck Mk Mk Mk Mk kx kv kx ko k ko k ko k kk kk s Sree Enae e rou KKK KKK KKK KKK KK KKK KKK KKK KK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK Orbit_azimuth_offset master slave LH ines Orbit_range_offset master slave 18 pixels timated_azimuth_offset master slave 193 787 lines timated_range_offset master slave 17 5939 pixels timated_azimuth_timing_error_lines master slave 3 lines timated_range_timing_error_pixels master slave 0 pixels timated_azimuth_timing_error_sec master slave 0 00178588 sec timated_range_timing_error_sec master slave 0 sec kk ck ck ck ck ck ck ck KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK ck ck kk ko kk kk ko kk Ck Mk ck kk ko ko kc kok ko nn n n o n End timing error NORMAL ck ck ck ck ck ck ck ck kk ck kk ck ck ck ck kk KK KKK ck ck ck KK KKK KKK ck kk kk ck ARA RARA k kk kk Mk ck Ck ck ko ko kc ko k ko Se Se AO dub dest dest usb det Teal In the logfile some additional informa
57. 2ell orbitbk cc int32 Ip2ell 115 compbaseline orbitbk cc void compbaseline 116 dumporbit orbitbk cc void orbit dumporbit 117 splineinterpol orbitbk cc matrix lt real8 gt splineinterpol 118 showdata orbitbk cc void orbit showdata 119 main orbitbk cc int82 main 120 eq1 doppler orbitbk h inline real8 eq1 doppler cn velocity cn dsat_P 121 eq2_range orbitbk h inline real8 eq2 range cn dsat P real8 rangetime 122 eq3 ellipsoid orbitbk h inline real8 eq3_ellipsoid cn P real8 semimajora real8 semiminorb 123 eq1 doppler dt orbitbk h inline real8 eq1 doppler dt cn dsat P cn velocity cn accerelation 124 shownumberofpoints orbitbk h int32 shownumberofpoints return numberofpoints 125 main processor c int32 main 126 handleinput processor c void handleinput int argc char argv input gen amp input general 127 usage processor c void usage char programname 128 fillproductinfo productinfo cc void productinfo fillproductinfo 129 readphase productinfo cc matrix lt real4 gt productinfo readphase 130 productinfo productinfo h productinfo multilookL 1 multilookP 1 rest 131 showdata productinfo h inline void showdata const show content 132 compinterfero products c void compinterfero 133 compcoherence products c void compcoherence 134 subtrrefpha products c void subtrrefpha 135 subtrrefpha products c void subtrrefpha 136 subtrrefdem products c void subtrrefdem 137 dinsar produc
58. 4298 SES 0 47 db s Sr OOD ASS SAS 8 268 2744 240 69 S Oo Wal 0 SaL 0 0 325520 19 65 SKIP SKIP ES E A MEAN NND Sa SS 0 41 0 14 0 18 150 86 184 45 MIS MAA ASIS AO eX De 52 DOS 0 07 GOL SU 12 21 494 14974 1938 241 00 225120516 0 42 0 2 p 0 14 280 63 147 84 20 3 Implementation The observation equations are given by the polynomial model y A x Q1 0 p 0 Ui 1h p Q pf 010 Ya 1 la po 3 din pj O01 MIN d 20 2 UN 1 iw pn Uy ce py aoa Where y contains the observed offsets in a certain direction denotes the location line number of the observed offsets in a certain direction p denotes the location pixel number of the observed offsets in a certain direction Qip denotes the unknown coefficients of the polynomial The data is rescaled to the interval 2 2 see Annex D so the normalmatrix is rescaled otherwise there could occur very high values for e g 25000 The least squares parameter solution is given by 59 APO y A 0 Ar Ne 20 3 Where q is the diagonal covariance matrix of the observations this matrix can be equal to identity or to the correlation values in version 1 CPM_WEIGHT card The coefficients are estimated by factorization of the matrix N The inverse of matrix N is also computed to check the solution stability and to compute some statistics A check number is given max abs N N I that gives a hint on the stability of the solution 60 Chapter 21
59. A kk kk ko Sk kv ko kv kv KKK Include flatearth No DEM source file d2 doris test final wanaF2835 dem Min of input DEM 92 Max of input DEM 17 Sal Data_output_file 42408_22735 demphase Data_output_format real4 First line w r t original master 3060 ast line w r t original master 8052 First pixel w r t original master TS Last_pixel w r t original_master 2710 Multilookfactor azimuth direction iL Multilookfactor_range_direction TL Number of lines multilooked 4993 Number of pixels multilooked 992 ck ck ck ck ck ck ck ck ck ck cock ck ck ck ck ck KARA ck ck ck ck ckock ck cock ck ck ck ckock ck ck ok ckock kk ck ck ck ock ck ko ck ck ck ck Ck ko ck ko ko kv Mk kx ck A AX x End comp refdem NORMAL ck ck 0k ck ck Ck ck ck ck ck 0k ck ck Ck ck kk ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck kk ck ck ck ck ck ck ck ck AAA AA AAA kk Sk Sk kv Sk kv kv kx ko The output in the logfile is more verbose specifying the results of the intermediate steps Also go over the standard out in case of problems with the SCREEN set to DEBUG level Figure 26 1 shows an example of a real valued phase map interferogram for a radarcoded DEM the Data output file This applies to the same scene as described in Section 23 0 4 26 3 Implementation The DEM reference phase is computed in two steps First the DEM is radarcoded to the coordinate systems of the master image Per DEM point the master coordinat
60. BAMLER NONE LINEAR QUADRATIC Experimental card Weight estimated offsets observations based on correlation in least squares solution weighting option Bamler was added in v3 16 and made the default recommended The theoretical precision of shift estimation using coherent patches is the basis of this weighting option CPM_MAXITER 10 Number of outlier to remove automatically based on outlier test The least squares adjustment is repeated until all tests are accepted or the max number of iterations is reached CPM_K_ALPHA 1 97 Critical value of outlier detection A higher value accepts more outliers This value can be found as the sqrt of normal distribution if you want a level of significance for the outlier test of 0 05 then look the value up under a half sided test Example input cards for this step e e comment COMPUTE COREGISTRATION PARAMETERS s CPM THRESHOLD 0 4 CPM DEGREE 2 CPM_WEIGHT linear none c CPM_WEIGHT quadratic none CPM_MAXITER 20 CEM _ PLOT NOBG 20 2 Output Description The plots are made if CPM_NOPLOT is not set The plotcpm uses a file CPM_DATA which is created in the working directory containing the data to be plotted The process control flag at the start of the products result file is switched to 1 at successful exit comp_coreg 1 Example of output of this step products result file KAEKKKKEKKEKKKKKKKKKKKKKEKKKKKKKKKKKEKKKKKKKKKKKEKKKKK
61. Bj l 1 B x Bh E v B 0 N Bj N Figure D 2 Definition of the baseline parameters a parallel perpendicular b horizontal vertical c length orientation Position 1 is the reference position Bj gt 0 when R gt Fi where R is the corresponding slant range The angle a is defined counter clockwise from the reference satellite 1 starting from the horizontal at the side of the look direction 146 denotes a pointwise multiplication lis the complex interferogram M is the complex master image S is the complex resampled slave image Ris the complex amplitude 1 reference phase The phase image of complex interferogram minus reference phase is defined as arctana Jimag treat D 16 Where arctan is the four quadrant arc tangent o is the phase image lis the complex interferogram Which is equal to with an ambiguity of 277 r M os OR D 17 The reference phase is defined as where r denotes the range from the master satellite to a point on the reference surface Ar An Which is the same as R M S D 19 Where M denotes the phase of a point situated on the reference surface of course in this definition the phase of the interferogram equals zero if there actually is no topography and M M The values of the real valued reference phase are stored in a 2d polynomial of certain degree The subtraction of the reference p
62. CK kk kk Sk Sk Sk Sk Sk Sk Sk Sk Sk Sk Sk Sk Sk kk Sk kk kk kk kk kk Sk kk Sk Sk kk kk Sk Sk kk Ck Ck Sk Ck AA KAR 118 33 3 Implementation Known are the heights of each pixel in the master line pixel system The point P x y z corresponding to a line pixel is computed with the 3 equations see Annex D in such a way that it lies on an ellipsoid of height h above the refernce ellipsoid When these coordinates are known the equations of Bowring are used to transform them to an ellipsoid system 6 A h The semimajor axis is denoted by a and the semiminor axis is denoted by b The squared first eccentricity by 2 72 2 at b The squared second eccentricity by e 1 e sin3 cos3 gt a p b2 arctang z a r b sin v cos v arctang z e b sin3 r e a cos3 arctan y x a v1 e sin m cos 119 83 1 Bibliography Arikan et al 2008 Arikan M van Leijen F Guang L and Hanssen R 2008 Improved image align ment under the influence of elevation In Fifth International Workshop on ERS Envisat SAR Interferometry FRINGEO7 Frascati Italy 26 Nov 30 Nov 2007 page 4 pp Bahr and V gtle 1991 Bahr H P and V gtle T 1991 Digitale Bildverarbeitung Anwendung in Pho togrammetrie Kartographie und Fernerkundung Wichmann Verlag Karlsruhe Curlander and McDonough 1991 Curlander J C and McDonough R N 1991 Synthet
63. Doris Delft Object oriented Radar Interferometric Software User s manual and technical documentation Version v4 02 Delft Institute of Earth Observation and Space Systems DEOS Delft University of Technology G TUDelft Preface This document describes the Doris Software for Interferometric SAR processing It is compliant with Doris v4 02 This manual contains technical documentation and information that is required for running the soft ware We try to be as complete as possible but some chapters may be very brief in the description Please report any incompleteness in the documentation or the source code The latest information on the Doris soft ware can always be found on the internet http enterprise lr tudelft nl doris Doris is freely available to the scientific community The conditions of use for the Doris software are as follows 1 Doris is a scientific purpose software and cannot be commercialized nor can parts or products of it be commercialized Parties interested in using Doris or its products for any commercial purposes are requested to contact Prof Dr Ramon Hanssen of DEOS r f hanssen tudelft nl 2 Our version of the software is the only official one Please do not distribute the Doris software to third parties instead refer to the Doris home page This in order to guarantee uniformity in the distribution of updates and information 3 Delft University of Technology is not responsible for damage of any kind cau
64. EER PPP PRP RP RPP RP PRP PP PP PRERREPRPRRPRPRRRRRPR PP PPP PPP RPP RPP RPP Pp PRPPPRPRPRPPRPRP PPP EE PRPPPRPRPRPPRPRP RPP POP The idea is to have a template directory with ascii kernels in it that then can be used in doris do not have much experience with this One could generate the different filters with Matlab For spectral method one may want to filter with an Hamming window If you have matlab paste the following to your terminal uses Matlab to generate the ascii file just an example please experiment yourself Moreover if you want to offer your ascii kernels for standard distribution of doris please email it to us matlab lt lt __EOFHD gt dev null Sum 395 filterfile filter hamming f standing hamming SIZE x ones 1 SIZE ones SIZE 1 xlying hamming SIZE fid fopen filterfile w torier nC ub Si Lo Na y SAD SES p for ii 1 SIZE noreste de tall Ga wie Y rode abate 2 E orcas ae atolo Nie y p end exit EOFHD And in Doris use the cards using filter in spectral domain in this case PF METHOD spectral PF IN KERNEL2D filter hamming PF_BLOCKSIZE 32 PF_OVERLAP 4 29 2 Output Description The process control flag for this step is switched on in the products result file 94 Figure 29 1 Magnitude of unfiltered complex in Figure 29 2 Phase of unfiltered complex interfer terferogram ogr
65. KKKKKKKKKKKKKKKKKKK _Start_coregpm KEKKKKKKKKKKKKKKKKKKKKKKKKKKEKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKK Degree_cpm il Estimated_coefficientsL 2 41088165e 02 0 0 1 48768713e 05 1 0 Sd se 05 0 il Estimated coefficientsP 58 3 11544442e 00 0 WoSVPSLSLOIS 06 il JL SHLILG1L2 05 amp 02 9 1 KKKKKKKKKKKKKKKKKKKKKKKK KK KK KK KK KK ck ckckckckckckckckckckckckckckckckck ck ck ck ck ck ck ck ckckckck k ck kk 0 0 End coregpm NORMAL Kckck ck ck ck ck ckckckckckckckckckckckckckckckck ckckckckckckckckckckck ck ckckckckckckckckckckck ck ckckckckck ck ck ck ck ck ck ckckckckck k ck kk In the logfile some additional statistical information is written The standard deviation of the estimates and the residuals after the least squares adjustment An ascii file CPM Data is created with some information for the plotcpm An example is shown below File CPM Data This file contains information on the least squares estimation of the coregistration parameters This info is used in the plotting scripts There are 10 columns containing Window number position L position P offsetL observation offsetP observation correlation estimated errorL errorP w test statistics for L P win posL posP Qi 11h OP Core eL eP wtstL wtstP 0 268 30 241 s O hq ALY 0 42 0 08 TLI 81 25 200 79 il 268 SES 241 00 o Si 0 35 0 01 Qs d 15 13 249 9 3 268 1048 241 00 3 38 0 46 0 01 O15 11 35 3790522 5 29 Ji 2 2
66. L i 1 i rm LH rit onli 1 HT I I eit ALII i I iH HATTE EE pP NO APER LAETI TT Pon om et TT THU TERT lt x zQ _ gt _ M w Figure 23 1 Phase image of complex interferogram The flat earth is clearly visible as a trend over the interferogram 23 0 4 NEW Note that for the optimal results i e to avoid an aliasing of spectras the images first have to be oversampled by a factor two before multiplication for an optimal result See sections on OVERSAMPLE card For a more modular approach in the new method it is not advisable to subtract the reference phase in this step If you do want to subtract the reference phase here then make sure you first run Doris to run comprefpha and then make a second run for step interfero After generation of the complex interferogram the reference phase can be computed by the new module comprefpha and subtracted by the new module subtrrefpha Also a reference height model can be computed and subtracted in future modules Figure 23 1 shows an example of a complex interferogram Only the phase is shown here We processed orbits 21066 and 1393 of frame 2781 Italy acquired at 26th and 27th July 1995 respectively ERS1 2 Tandem mission The parallel baseline is about 35 meters which implies a height ambiguity of about 270 meters Clearly a large trend caused by the flat earth is present
67. Makefile Follow the directions on the screen Compile Doris and install the executables in directory Source_new inspect edit the Makefile make this compiles the code make test this should give the version number of Doris make install uses usr local bin by default Also the bin utilities are installed Make sure the Installation directory is in your path For t csh users it should be in your t cshrc file startup file Add it with a csh command like set path usr local bin path B 1 2 Installation of the SARtools We also need to compile the SARtools and ENVISAT TOOLS programs There are Makefiles in the sub directories provided with default installation directories If you selected another installation directory than usr local bin please change that in the two Makefiles What you have to do to install these utilities is O oc A C Iv cd SARtools review the Makefile make n check what happens make compile software make n install check if this is what you want make install B 1 3 Installation of the ENVISAT tools And for the ENVISAT TOOLS o a RK wo nm cd ENVISAT_TOOLS review the Makefile make n check what happens make compile software make n install check if this is what you want make install We did not automate this because of the complexity while simply editing two Makefiles should not be a big problem 125 B 1 4 I
68. P OUT slave raw Filename of the raw data output file S DBOW linelow linehi pixellow pixelhi Slave output window You can make a cutout of the image with this card If card omitted it defaults to total image line pixel 1 refers to the the first line pixel S DBOW GEO lat 0 lon_0 height width Slave output window Alternative to and overrides normal DBOW card You can make a cutout of the image with this card latitude of the center pixel of the desired crop longitude in decimal degrees WGS84 system of orbit then height width in pixels For approximately square areas heights should be a factor 5 of width for ERS 32 Chapter 12 S OVS In this chapter the processing of step S OVS is described It is the same as step M OVS but then for the slave image See chapter 8 M OVS for more information on this step 12 1 Input Cards S OVS OUT Slave ovs raw Filename of the oversampled data S OVS OUT FORMAT ci Output file format S OVS FACT RNG 1 oversampling factor of output image in range pixels S OVS FACT AZI 1 oversampling factor of output image in azimuth lines S OVS KERNELSIZE 16 Kernel size sinc used for oversampling 33 Chapter 13 COARSEORB This chapter describes the processing step COARSEORB In this step the coregistration based on the orbits of slave and master is computed with an accuracy of about 30 pixels precise orbits This is a fast way to get the coarse offsets Before the FINE coregistration
69. Possibly in your configuration there is a file complete tcsh that is sourced from the cshrc If you are not using tcsh you cannot use them Find out by the commands who am i and finger complete doris c c h q ver A n h lt search term gt n f in IN doris x complete cpx2ps c wfqes lLELpPTFezoGCgKSUVZhm 7 n f ci2 cr4 cr8 r4 n q normal mag phase real imag Y n T title gt n c cool copper gebco gray haxby hot jet no green polar rainbow red2green relief topo sealand split wysiwyg n m X XY Y complete cpxfiddle c w f qoeslLpPSMmcVh N n c filename gray jet hot cool bert n f ccl cucd ci2 ci4 cr4 cr8 K n q normal mag phase real imag n m X XY Y n o ascii float sunraster uchar 143 Annex D Definitions In this annex a number of definitions as used by Doris are described In Section D 2 the baseline repre sentations are described while in Section D 3 the definition of the interferogram is described Section D 4 describes the definition of the polynomials And in Section D 5 the use of the pulse repetition frequency and range sampling rate are discussed Section D 6 describes a system of 3 equations which is frequently used to compute the position on the ellipsoid for a certain line pixel The way the orbits are interpolated is described in Section D 7 Finally sectio
70. RESAMPLE In this chapter the resampling or interpolation of the slave image on the master grid is described The slave image is resampled reconstruction of original signal from the samples by correlation with interpolation kernels in space domain accordingly to the transformation model from the step COREG_PM and optionally DEMASSIST This model states with sub pixel accuracy which points of the slave correspond to the master grid Note This step may be fairly time consuming The spectrum in azimuth can be centered at zero before resampling and shifted back to its original Doppler centroid frequency afterwards This is required since the spectrum of the kernel function is centered at zero See also Geudtner 1996 This shifting has been implemented in release 3 0 but not in prior versions The polynomial described by the strings in the slave result file Xtrack f DC constant etc is used If the spectrum should be shifted use the card RS SHIFTAZI New in v3 4 is that the azimuth kernel is shifted to the Doppler centroid not as before the dataspectrum to zero and back This is made default To assess the quality of the resampling the resampled slave image can again be coregistered step FINE onto the master This should yield offset vectors that are normally distributed with zero mean The slave image could also be resampled in steps first resampling it by a first degree model then again with a higher degree model 21 1 Input C
71. SS SS You should add the Bin directory to your path The utilities in the Bin directory are there for your convenience and you may edit them the way you prefer Possibly you need to make these files executable by the command chmod 755 Bin x Bin run Generate input and shell for running Doris Bin helpdoris Summary of input keywords print with helpdoris p Bin plotcpm Step coregpm GMT show estimated error offset model histograms etc Bin plotoffsets Step fine GMT to show offset vectors thresholded on correlation Bin viewanddel Step coregpm background call to gv thereafter delete dummy file B 7 List of routines description The information in this section may be out of date The list is generated with the ctags command bin ctags u x x ch grep v matrixspec grep v matrixbk grep func cut c1 20 39 600 sqr constants h inline int16 sqr int16 x return x x sqr constants h inline int32 sqr int32 x return D sqr constants h inline uint sqr uint x M Xx X sqr constants h inline real4 sqr real4 x return x x sqr constants h inline real8 sqr real8 x return x x in constants h inline real8 in cn P const scalar product r P in Q out constants h inline cn out cn P const cross product cn r P out Q dist constants h inline real8 dist cn P const distance d P dist Q min constants h inline cn min cn P const cn rzP min Q norm2 constants h inline re
72. Sk S ke ke Sk kc ko k ck ck ckckok We noticed that if the precise orbits are not long enough not enough time before and after first last line this results in a wrong reference phase for obvious reasons Interpolation near the end of the data points is not very good with cubic splines This can be solved by using more orbital data points after the last line of the Scene see cards for step M PORBITS Figure 25 1 shows the result of subtracting the reference phase polynomial from the interferogram Fig ure 23 1 The same scene as described in section 23 0 4 again multilooked now by factors 4 and 4 resulting in total multilooking of 40 and 8 which agrees on the terrain with a resolution of about 160 meters square This step can be mis used to correct for residual orbital fringes if you know what you are doing To do this first count the number of fringes you want to remove from the interferogram Then edit the products result file and create a section for the step COMPREFPHA In this output section simply define a polynomial that describes for example a linear trend in range of say 2 5 fringes Then run this step and doris will not know that it is not the reference phase polynomial that is subtracted from the interferogram but an additional correction polynomial that has been inserted by hand 81 Figure 25 1 Phase image of complex interferogram The flat earth is subtracted leaving dominantly topo graphic fringes
73. TALLDIR make this should compile the code make install this installs in usr local bin 2 Make sure the Installation directory is in your path For tcsh users it should be in your cshrc file startup file Add it with a csh command like set path usr local bin path C 2 2 run script This script can generate template input files and directories and therefor generally speeds up the processing If it is installed in a Bin directory in your path it can be run from several project directories for uniform processing It uses the environment variables EDITOR and PAGER if set The basic idea is to start with run g to generate the input and then to to edit default with editor vi or EDITOR the first input file run el After you saved the file located in directory Inputfiles type run s1 to process call to doris software the first input file The output that goes to stdout can be viewed with run v1 It is redirected to a file in directoy Outinfo To view the result files that are created by Doris use run r1 master run r2 slave run r3 products run r4 logfile This should be repeated for other steps Of course this run file is only a helping hand not the solution to all your problems Be careful Template values aren t always the best settings Feel free to improve it Usage run s e v step f inputfile r file id d g M master S slave
74. TPO bl 4B2 Q 47 Lr 60 lt E6 1499 a eee RT LN e 1499 ep E 6 289 B1 42 H3 4o 5 n6 us ez Ba p en a Sp B1 T T T T T T T T T 399 798 1197 1596 1995 2394 2793 3192 3591 3990 GUT EEE Range Figure 20 4 Plot produced by plotcpm for the second run The estimated offsets are plotted here normal ized together with a 90 degrees rotated w test as ellipses 55 Azimuth_ direction Correlation 0 0 0 1 0 2 0 3 0 4 0 5 0 6 0 7 0 8 0 9 1 0 0 0 0 1 0 2 0 3 0 4 0 5 0 6 0 7 08 0 9 10 Correlation d MIT Feb 16 11 51 35 2000 Figure 20 5 Plot produced by plotcpm for the second run The absolute error estimated offsets minus observed offsets are plotted for azimuth direction Range_direction Correlation 0 0 0 1 0 2 0 3 0 4 0 5 0 6 0 7 0 8 0 9 1 0 0 5 0 5 0 4 0 4 NS Dos o3 L o 4115 o 5 0 2 eo a 0 2 E 0 1 m E P MN 0 1 0 0 f3 0 0 00 01 02 03 04 05 06 07 08 09 1 0 EM resto 151362000 Correlation Figure 20 6 Plot produced by plotcpm for the second run The absolute error estimated offsets minus observed offsets are plotted for the range direction 56 Location _windows wtests Fanga mg m wp a m Figure 20 7 Plot produced by plotcpm The magnitude is plotted in the background d i Tey y ase yl 20 1 i 0 j 0 Perhaps one might do the resampling with a lower quality polynomial and thereafter do the fine coregistration initial offsets 0 0
75. To Rk ee eR eoe Ke Re KR ec X ko ko KR eo ce ee hk ae ek e ee Kok kk ek ete He ak ie dee ak ak Fe He eK Start subtrrefdem okckckckockckckckockckck kockockck ckokokckck o kock EXA kokcokck kockckok ck RARA kockockok k ck ck ckck kck kck AEREA RARA RAROS Additional azimuth shift T Additional_range_shift 2 87 Figure 27 1 Phase image of complex interferogram The phase of the reference DEM is subtracted flat earth was already subtracted in subtrrefpha leaving residual topographic atmospheric and errorous fringes Data_output_file Outdata cint minrefdem raw Data_output_format complex_real4 First_line w r t original_master DAL Last line w r t original master 14964 Es E O nat os en in Last pixel w r t original master 3998 Multilookfactor azimuth direction 40 Multilookfactor range direction 8 Number of lines multilooked 368 Number of pixels multilooked 499 ckCkckckckckckckckckckckckckckckckckckckckckckckckckckckckckckckckckckckckckckckckockckockckok ckock ck ckokckckckckckckckck ck ck ck ckck kc kk x End subtrrefpha NORMAL ck ck kc ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck Ck ck ck ck ck kk ck ck ck ck ck ck KK ck ck ck ck ck ck ck ck ok ck ck RARA ok kk ck ko Sk ck Sk kv kx ko ko ko Figure 27 1 shows the result of subtracting the reference phase of the radarcoded DEM from the interferogram See Figure 25 1 It can be seen the number of topographic fringes is reduced though there still
76. You can create such a file with appropriate dimensions using Matlab for example Alternatively the much faster Unix way would be along these lines First compute the height number of lines of the dummy file echo 2105 10013 1 100 be i 110 5 then the width number of pixels echo TO00 301 1 2 be i 100 ie the file should be 110 lines by 100 pixels of 4 byte Now create the file using dd dd if dev zero of Outdata dummy height raw count 110 bs 400 33 2 Output Description The process control flag is switched at successful exit geocoding il An example of the output in the products result file KKKKKKKKKKKKKKK KK KK KK KK KK KKK KKK KKK KKK KKK KKK KKK KKK KKKKKKKKKKKKKKKKKEK Start geocode kckCckCckckck ck ckckckckckckck ck ckckckckckckckckckckckckckckckockckckckckckckckckckckckckckckckckckckckckck ck ck ck ck ck ck ck ckck ck ck k ck kk Data output file hei slant2h Outdata hei ambi Data output file phi Outdata phi raw Data output file lamda Outdata lambda raw Data output format real4 First line w r t original master TOOL Last_line w r t original_master 205 kims t jambrel Most onigi nalimas EEn 501 Last pixel w r t original master 700 Multilookfactor azimuth direction AO Multilookfactor_range_direction 2 Re Ne Fe eoe eode ox oe oe Pee te te Ke ae Te ke e te te He hee te Xe she o eoe Ie at ah te te Ute oe ee Ie ete ee Ie aie ee Te a ae te ie Hee x End geocode NORMAL kk
77. al master 2713 Number of lines 5000 Number of pixels 1000 Deltaline slave00 dem 194 341 Deltapixel slave00 dem 17 7004 51 Deltaline slaveON dem 194 082 Deltapixel slaveON dem 18 8867 Deltaline_slaveN0_dem 194 479 Deltapixel_slaveN0_dem 17 8064 Deltaline_slaveNN_dem 194 232 Deltapixel slaveNN dem 18 9409 CK CK CC Ck CK CI CI CK C Ck CK CC CC Ck Ck Ck CI CIC Ck Ck CK CI CK CI C CC Ck Ck C Ck C CK CI C I I E E Kk ko Sk Sk Kk kx kx x M x End dem assist NORMAL ACkCkCck ckCckck ckckck ck ck ck kck ck kck ck ck ck ck k ck ck kk ck k ck k ck ck ck ck ck ck k ck ck ck ck ck ck ck ck ck ck ck ck kk ck ck k ck k ck ck k k k kk kk The output in the logfile is more verbose specifying the results of the intermediate steps Also go over the standard out in case of problems with the SCREEN set to DEBUG level 19 3 Implementation The DEM assisted coregistration is estimated in two steps First the DEM is radarcoded to the coordinate systems of the master and slave image Per DEM point the master coordinate and the offset between master and slave is saved to a file Both the master coordinates and offsets are real valued Second the offsets are interpolated to the integer grid of master coordinates A linear interpolation based on a Delaunay triangulation is used The software package Triangle for the Delaunay triangulation is kindly made available by Jonathan Shewchuk Shewchuk 1996 Shewchu
78. al8 norm2 const n P norm2 norm constants h inline real8 norm const n P norm normalize constants h inline cn normalize const cn R P normalize angle constants h inline real8 angle cn A const angle A angle B 0 pi 134 eccist sqr constants h ecc2nd_sqr constants h disp constants h lines constants h pixels constants h operator constants h operator constants h operator constants h pol2xyz conversion c xyz2pol conversion c xyz2ell conversion c xyz2ell conversion c ell2xyz conversion c deg2rad conversion h deg2rad conversion h rad2deg conversion h rad2deg conversion h line2ta conversion h pix2tr conversion h pix2range conversion h ta2line conversion h tr2pix conversion h cr4toci2 conversion h coarseporbit coregistration coarsecorrel coregistration coarsecorrelfft coregistration corrfft coregistration distributepoints coregistration getoffset coregistration finecoreg coregistration coherencefft coregistration coherencespace coregistration coregpm coregistration getofffile coregistration cc4 coregistration cc6 coregistration ts6 coregistration ts8 coregistration ts16 coregistration rect coregistration tri coregistration resample coregistration rangefilter filtering c rfilterblock filtering c phasefilter filtering c goldstein filtering c smooth filtering c smooth filtering c spatialphasefilt filtering c convbuffer filtering c phasefilterspectral filtering c spectr
79. alfilt filtering c azimuthfilter filtering c blockazifilt filtering c slant2hschwabisch geocode c slant2hambiguity geocode c slant2hrodriguez geocode c geocode geocode c printcpu ioroutines c inittest ioroutines c doinitwrite ioroutines c initwrite ioroutines c updatefile ioroutines c getanswer ioroutines c readres ioroutines c updateprocesscontrolioroutines c checkprocessing ioroutines c checkrequest ioroutines c fillcheckprocess ioroutines c fillprocessed ioroutines c filelines ioroutines c lines inline void ecc ist sqr first ecc inline void ecc2nd_sqr second ecc inline void disp const show content inline uint lines const return number of inline uint pixels const return number of pixels window amp operator Operator operator pol2xyz xyz2pol xyz2ell xyz2ell ell2xyz bool bool void void void void void inline inline inline inline inline inline inline inline inline real8 real4 real8 real4 real8 real8 real8 real8 real8 window X window X const window X const real8 real4 deg2rad deg2rad rad2deg real8 return return 2 Pll 7 180 3 s Pl y 180 7 return e 180 s rad2deg real4 return s 180 s line2ta real8 line real8 tal real8 prf pix2tr real8 pixel real8 tr1 real8 rangesamplingratex2 pix2range real8 pixel real8 tr1 real8 rangesamplingratex2 ta2line real8 azitime real8 tal real8 prf XXX x p
80. am filtphase in And in the same result file a section will be added like except if PF IN FILE is specified Xe cox ee Ke seo ee ee Fee oe oe ok ARA to ie Sec Koo oko eS le ee ic ee a Ie Xo e XS ox Start filtphase Ck ck ck ck ck ck Ck Sk ck Ck ck kk Ck Sk ck Ck ck kk ck kk ck kk Ck kk ck Sk ck kk ck ck KKK ck kk kk ck kk kk ck kk Sk Sk Sk kk KKK KK KK Method goldstein size alpha overlap 32 0 5 4 Input_file Outdata cint srp raw Data_output_file bmi a Oo BENE Data output format complex real4 First line w r t original master 1073 Last line w r t original master a 3012 First pixel w r t original master 148 Last pixel w r t original master 985 Multilookfactor azimuth direction 10 Multilookfactor range direction 2 Number of lines multilooked S23 Number of pixels multilooked als KKKKKKKKKKKKKKKKKKKKKKKKKKAKXKKKKKKKKA KKK KKK KKK KKK KK KK KA KK ko ko ko ko ko ko End_filtphase _NORMAL Ck ck ck ck ck ck KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK ko kk Sk kv Sk kk Sk kc k ck ck ko 29 3 Implementation 29 3 4 spatialconv The complex interferogram is convoluted with a kernel by FFT s The card PF_KERNEL specifies the 1D kernel The 2D kernel is computed as PF KERNEL PF KERNEL e g for a 3 point moving average 1D 95 Figure 29 3 Magnitude of filtered complex in terferogram Method spatialconv A spatial convolution with a kernel 1 4 9 4 1 was used Clear
81. aneously for all alphas because the normal matrix and factorization remains the same but somehow our cholesky routine introduced an error which is probably caused by using fortran in c so we just solve 3 seperate times with the same factored normalmatrix 4 Evaluate for all points at line pixel with ok unwrapped phase the 2d polynomial to obtain the coefficients of the height phi function Then evaluate the 1d polynomial to obtain the height For all l p with ok unwrapped phase e Compute the alphas For betas appropriate to alphaO d D w M Boss ip 32 28 i 0 j 0 e For betas appropriate to alpha1 d D a x fij dpi 32 29 i 0 j 0 e repeat computing alphas until you have them all 1dD 1 e Compute the height 1dD h Y ad 32 30 i 0 32 4 Comparison of the methods Here a simple test is described that was performed to see the differences between the methods A unwrapped interferogram was obtained of the Veluwe Holland by processing the bottom half of the Tandem images 3512 ERS2 and 23185 ERS1 The interferogram was multilooked by 40x8 resulting in 367 lines and 610 pixels There were about 3 fringes Baseline B 185m a 3 Bj 62 B 173 Bj 184 B 2 The processing was done with a debugger version of the Doris software so the cpu times are not really representative for the performance of Doris table 32 1 shows some processing parameters Figure 32 2 and 32 3 show plots for these thr
82. ards RS_METHOD RECT TRI CC4P CC6P TS6P TS8P TS16P KNAB6 KNAB8 KNAB10 KNAB16 RC6P RC12P Select kernel for interpolation a simple step function nearest neighbor or a linear interpolation tri or cubic convolution kernel 4 or 6 point or a knab sampling window using a default data oversampling factor or a Raised Cosine best kernel 6 or 12 points or a truncated sinc 6 8 or 16 point can be used RS OUT FILE s resampled raw Output filename of resampled slave cannot be equal to the input file name RS OUT FORMAT CR4 CI2 Output format of resampled slave complex_real4 or complex_short Same as SLC input format this causes an error of maximum about 1 percent 61 RS DBOW linelo linehi pixello pixelhi Data base output window in master coordinate system for slave to make a cutout If card is omitted it defaults to the overlap between master and slave and corrected for half the kernel size where no interpolation is possible For stacking of interferograms on top of each other use the coordinates of the master after cropping In this way all interferograms are automatically aligned If the slave image is smaller than the window the pixels are set to 0 RS_DBOW_GEO lat_0 lon_0 height width Output window specified in latitude longitude decimal degrees WGS84 Alternative to and overrides normal RS_DBOW card The latitude and longitude of the central pixel of the desired crop are specified toget
83. aseline is approximately 30 meters And corrected for the phase of the reference body o fa Bjo 31 8 For the defo pair 1 3 denoted with a prime similar equations follow Deformation in the line of sight range that occurred in between the acquisitions is denoted by Ar Ar os 31 9 4m A positive Ar implies deformation in the B direction away from the sensor i e subsidence The phase of this interferogram is EU d r ra A V Ar 31 10 Combining the expressions for the interferometric phase for the topo pair 31 7 and defo pair 31 10 yields B 4m Il o Pp A 31 11 The problem here is that the true parallel baselines are unknown The actually wrapped phase of the deformation interferogram corrected for reference phase is defined as 4 IB Bio Av 31 12 z Fp sin 0 o B sin 0 o Ar I B sin 00 B sin f Ar 105 where 3 4 a Using the approximation for small 6 which is about 1 or 0 0175 rad for terrain height differences of 5 km sin B 60 sin B cos 00 cos f sin 00 sin 8 60 cos B 31 13 it follows from equation 31 13 that the flat earth corrected phase equals 4 Z B sin 6 66 cos f B sin 8 Ar 81 14 4 ep cos 3 Ar An Am 66B 9 A xn Ee Gad The corrected phase for the topo pair equals 50B o using the same approximation
84. ave been done for the slave image however we have chosen to use the master slave fine coregistration because this is probably more precise for the relative timing Therefore an error in the master timing will propagate to the slave timing however the relative DEM position is consistent 18 1 Input Cards RTE THRESHOLD 0 4 Threshold for correlation value to use estimated offset of step FINE in estimation of the correlation based offset This depends on the size of the window during FINE Estimated coherence using small windows are more biased towards 1 0 so a higher threshold is better For window size 64 64 a threshold 0 2 seems OK The plotoffsets script can be used from the prompt to figure out a good threshold value RTE MAXITER 10000 Maximum number of iterations in least squares adjustment and testing procedure The least squares adjustment and testing is repeated until all tests are accepted or the max imum number of iterations is reached If the maximum number of iterations is reached before acceptance of all tests this is an indication that the fine coregistration failed or is of bad quality In this case either run the FINE step again obviously with different settings or increase the RTE K ALPHA card to accept a lower quality coregistration RTE K ALPHA 1 97 Critical value of outlier detection A higher value accepts more outliers This value can be found as the sqrt of the normal distribution Typical values are 1 97 and 3
85. ay triangulation is kindly made available by Jonathan Shewchuk Shewchuk 1996 Shewchuk 2002 See also http www cs cmu edu quake triangle html The interpolated topographic height and look angle are then used to compute the local incidence angle 6 for each point Using the local incidence angle we obtain the synthetic amplitude by synthetic amplitude sin 070c 1 6 1 The 1 is applied to obtain positive numbers Finally the obtained simulated amplitude per master pixel are saved to a file and used in the M_TIMING step to obtain absolute timing error of the master acquisition 24 Chapter 7 M_TIMING In this chapter the processing of step M_TIMING is described The absolute timing error between DEM and the master is computed using the simulated amplitude and the result of the coarse coregistration Therefore this step can only be run after the steps M CROP and M_SIMAMP The timing error in azimuth as well as in range direction is estimated During the coregistration the master is aligned with respect to the DEM based on the simulated amplitude image producing a single offset for the whole image Because the resolution of a DEM is typically coarser e g SRTM3 90 m than that of a radar image and the radar ground resolution differs between azimuth and range direction e g 4 m in azimuth and 20 m in range for ERS1 2 and Envisat the sensitivity in coregistration differs between the two directions measured in res
86. b package are used to obtain the points Note that getorb also interpo lates based on 30 second ephemerides We would like to test if setting the data interval to e g 30 second gives better results A test can be easily performed for the computation and modeling of the reference phase Assume that this phase can be accu rately modeled by a 2d polynomial of degree 5 Now first let the precise orbit be given with a data interval of 1 second Use step REFPHA to model the reference phase based on 501 points distributed over the scene Next let the precise orbit be given with a data interval of 30 seconds and again model the reference phase In the log file some statistics on the error of the model w r t the computed reference phase is given which can be used to find out which orbit gives a better model In both cases use at least 6 points before and 6 points after the last point in the frame 152 The interpolation is done as follows compare with numerical recipes in c splint routine First the piecewise polynomial coefficients are computed by solving a tridiagonal system and stored For interpolation the correct coefficients interval are read and the polynomial is evaluated Because we know that in our situation with getorb we always have ephemerides with a constant time interval we could speed up the computations Also the fact that this interval equals 1 can be easily exploited However we decided not to exploit these features because we l
87. between 90 90 and smallest longitude between 180 180 26 1 Input Cards CRD IN DEM filename filename of input DEM File is assumed to be stored in a raster Major row order from North to South line by line See also internet links at Doris home page for available DEMs CRD IN FORMAT 12 12 BIGENDIAN R4 R8 format of input DEM on file signed short for gtopo30 or real4 or real8 input ma trix is raw binary data w o header endianness of host platform is assumed except for I2 BIGENDIAN CRD IN SIZE 6000 4800 Number of rows and columns of input DEM file Default is set to tile w020n90 DEM CRD IN DELTA 0 00833333333333333333 deltalon 83 Grid spacing of input DEM in decimal degrees latitude longitude Default is equal gridspacing default set to tile w020n90 DEM CRD IN UL 89 995833333333333 19 995833333333333333333 Coordinates of UL upperleft corner in decimal degrees latitude 90 90 longitude 180 180 Default is set to tile wO20n90 DEM It is interpreted as max latitude min longitude in source CRD IN NODATA 9999 Identifier to ignore data in input DEM with this value Default is set to tile w020n90 DEM CRD INCLUDE FE OFF ON If this card is switched on the reference phase of the DEM is computed including the flat earth term Otherwize the phase is computed with respect to the ellipsoid yielding only topographic phase CRD OUT DEM filename Request optional debug output to float file
88. but also some topographic features can be seen In the frame the elavation ranges from zero to 1400 meters The original SLC images were cut out to 20000 lines by 4000 pixels The interferogram is multiooked by factors 10 in azimuth and 2 in range which yields a dimension of 1475 lines by 1997 pixels 23 1 Input Cards INT_OUT_CINT filename filename of output datafile for complex interferogram of step interferogram one of INT_OUT_ is mandatory INT OUT INT filename filename of output datafile for real interferogram of step interferogram one of INT OUT is mandatory INT MULTILOOK 51 74 multilookfactor if no multilooking is desired set this to 1 1 If the reference phase is not subtracted in this step be carefull not to multilook too much in this step In step subtrrefpha again a multilook card is present where one can multilook by factors 2 2 for example Example input section c comment _ product generation___ e ANO UN AN Output cint raw optional INTAOUT INT Output int raw optional CNO E Output flatearth raw optional INT_MULTILOOK i0 2 line pixel 23 2 Output Description At successful exit the process control flag is switched on interfero il The output looks like in the products result file eR eoe Ke He He Ke Ke ee se oe Fe ke ke He Kk ee A ie ee ake koc Hk ak ok dk ake oon te i ae le eee eR Xo ee XX es Hee ok x Start interfero kk Ck ck ck ck ck ck ck ck ck ck c
89. cols scalingfactor The next numlines lines contain the filter numlines and numcols must be odd centered for method spatialconv For method spectral they may be even the zero frequency is located at position kernelsize 2 1 starting at 0 The values of the kernel are multiplied by the scale factor The kernel is not normalized in any other way The output matrix has a zero valued edge of size floor kernel 2 Example of the cards for this step e e comment PRHASEFTLT_ a c PF_METHOD spectral c PF IN KERNEL2D proto myfilter CREPES HOCK SiH c PF_OVERLAP 4 E c PF_METHOD spatialconv c PF_KERNEL 5 3i 3L 3L TL db c PF IN KERNEL2D proto myfilter c PF_METHOD goldstein PF_IN_FILE Outdata cint srp raw 323 PF_ALPHA 0 5 PF_KERN apd beri eec PF_OVERLAP 4 EF UBLOCKSIZE 32 A simple example of PF IN KERNEL2D in an ascii input file use this example with method spatialconv S 5 0 05 93 Oe i dn 3 3 l 2 o db d ed 4 LEO A second example of PF IN KERNEL2D ascii input file One could use cards PF BLOCKSIZE 32 and PF OVERLAP 4 PF METHOD spatial i5 305 p pa f PPP PPP RPP PPP App PRPPPRPRPRPPRPRP RPP PP PRPPPRPRPRPPRPRP RPP EB PRPPPRPRPRPPRPRP RPP EE PPPPRPRPRPPRPRPRPRPRPE PO PRP PPP RP RPP PPP PP PPO PPPPRPRPRPPRPRPRPRPRPE RB PRPPPRPRPRPPRPRP PPP EB PRPPPRPRPRPPRPRP RPP PO PRPPPRPRPRPPRPRP RPP
90. d with the slave spectrum yielding loss of coherence in the interferogram 15 1 Input Cards AF BLOCKSIZE 1024 Length of fft per buffer in azimuth direction In general the larger the better AF_OVERLAP AF_BLOCKSIZE 8 Half of the overlap between consecutive bufferes in azimuth direction Partially the same data is used to estimate the spectrum which might have certain advantages However it has not been studied yet if taking an overlap is requird Setting this card to 0 is fastest AF HAMMING 0 75 The weighting of the spectrum in azimuth direction The filtered output spectrum is first de weighted with the specified hamming filter then re weighted with a newly centered one If this parameter is set to 1 no weighting is performed For more information see 22 3 3 AF OUT MASTER master afiltered Output file name for the master image 39 AF_OUT_SLAVE slave afiltered Output file name for the slave image AF_OUT_FORMAT cr4 Format of outut data Either complex real4 cr4 or complex shorts ci2 An example of the input file save general cards PROCESS motora PROCESS g irc E JL lf E comment _ AZIMUTH FILTERING m E iy c AF_METHOD AF_BLOCKSIZE 1024 fftlength each column AF_OVERLAP 64 ff i AF HAMMING 10 S AF OUT MASTER vieles Assia iNe AF OUT SLAVE Outdata 21066 azifilt AF OUT FORMAT cu 15 2 Output Description In the process control array the switch for azimuth filtering is
91. dependent of degree d d 0 Apo 1 d 1 Ajo Ao1 2 3 d 2 A290411 Apa 4 5 6 d 3 A30491 412403 7 8 9 10 Thus the number of coefficients unknowns in least squares estimation equals for a 2d polynomial of degree d 1 3 d 1 c d 1 D 27 And the degree of a polynomial with N coefficients is equal to d 5 int32 V 8N 1 1 D 28 D 4 1 Computation of coefficients Suppose we have 2d data f l p I 1 25000 p 1 5000 and we want to estimate a 2d polynomial D 26 with these data The system of equations looks like F l1 p1 1h p hm pi ES Flo p2 1 la po l lap pj en azo D 29 i H O11 f in pn 1 ly pw Inpn pr aod The convention used in the Doris software is that we first normalize the data to avoid numerical instabilities see source utilities hc The maximum coordinates are that of the original master stored in the result file of the master image typically 25000 lines and 5000 range pixels The coordinates are rescaled to the interval 2 2 l a la b nx i m D 30 Another way perhaps a better one would be to make the data zeromean unit standard deviation The estimated coefficients thus correspond to the normalized data For evaluation the data has to be normalized by the same factors a b Normally the information from the master originalwindow linelo etc are used e g for the coregistration and the reference surface polynomial These numbers can be found in the mast
92. dinates available you will have to create one This means that you will have to edit the products result file and create a SLANT2H section see Chapter 32 for a description of this section in the products result file The section will look something like 117 EXAXEXEXEXXENEK ERRE RARA RARE RARA EA AE AAA ERA AA ERE kck kockok k ko ota Ss lanos kk CK kk kk Sk kk Sk S KK KK KK KK KK KK KK KK KK KK KK KARA Sk Sk Sk kk kk kk kk Ck Ck Ck Ck Ck A KA A A AA KAR Method schwabisch Data_output_file Outdata dummy_height raw Data_output_format real4 First line w r t original master 1001 Last line w r t original master 2I First pixel w r t original master EDI tasr sles A o A dS EHE 700 ultilookfactor azimuth direction LO ultilookfactor_range_direction 2 Ellipsoid name a b WGS84 6 37814e 06 6 35675e 06 Ck Ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck Ck KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKKKKKKKK KKK KKK x End slant2h NORMAL ck ck CckCckckckckckckckckckckckckckck ck ckckckckckckckckckckckckckck ck ck ckckckck ckck kckckckckck ck ck ck ck ck ck ck ck ck k ck ck ckck ck ck k ck kk Also set the pcf to 1 on top of the products result file If you have an external DEM you can compute the required file Outdata dummy height raw using step COMPREFDEM The dimensions and multilooking can be copied from the interferogram section If your area is flat you may want to use a dummy file filled with zeros
93. dvolume void readleader void readnull void readdat void writeslc void unwraptreeframon void getorb void convertgetorbout void solve33 matrix lt real8 gt solve22 uint nextpow2 real8 polyval real8 polyval matrix lt real4 gt polyval real8 polyvalid void normalize void normalize void BBparBperp real8 8B void BBhBv int32 Btemp void BalphaBhBvBparBperpTheta real8 amp Bpar real8 amp Bperp inline bool iseven int16 w return w 1 2 inline bool iseven int32 w return w 1 2 inline bool iseven uint w return w 1 2 inline bool isodd int16 w return w 2 inline bool isodd int32 w return w 2 inline bool isodd uint w return w 2 inline bool ispower2 uint w inline int32 Ncoeffs int32 degree inline int32 degree int32 Ncoeffs inline real4 remainder real4 number real4 divisor inline real8 remainder real8 number real8 divisor inline real4 sinc real4 x inline real4 rect real4 x inline real4 tri real4 x inline real8 onedecimal real8 x inline real4 onedecimal real4 x inline matrix lt real4 gt myrect const matrix lt real4 gt amp X inline matrix lt real4 gt myhamming inline real4 interpbicubic inline real4 normalize real4 data real8 min real8 max inline real8 normalize real8 data real8 min real8 max 137 Annex C Utilities In this annex the additional software that aids running Doris is described This additional software is not required to
94. e C cptfile gKSUVZ h elp cpxfile For more info type cpx2ps h more 140 C 2 5 phasefilt doris Program to perform phase filtering from the prompt using Doris Several methods can be used For more info type phasefilt doris h PROGRAM phasefilt doris filter mph file using Doris SYNOPSIS phasefilt doris numlines o filtered a 0 25 s 2 e 3 b 32 2m goldstein k file infile C 2 6 flapjack Program to make integer linear combinations of interferograms PROGRAM flapjack pixelwise complex integer multiplication of a float complex file To be used to make linear combinations SEE ALSO cpxmult USAGE flapjack infile1 factor 2 EXAMPLE flapjack cint raw 3 C 2 7 cpxmult Program to subtract phase in 2 complex files PROGRAM cpxmult subtracts or adds phase of two complex float files USAGE cpxmult infile1 infile2 outfile add 0 infile 12 contain complex values a ib outfile contains a1 ib1 conj a2 ib2 If add 1 then phase is added not conj EXAMPLE cpxmult cint raw cint2 raw subtract raw C 2 8 cpxdiv Program to performs division on 2 complex files Phase is subtracted by default Program cpxdiv divides two given complex float files USAGE cpxdiv infile1 infile2 outfile cnj 0 infile 12 contain complex values a ib outfile contains a1 ib1 a2 ib2 If cnj 1 then phases are added outfile contains a1 ib1 co
95. e dl The output section in the result file will resemble the following kk kk kk Sk Sk Sk Sk Sk KK KK KK KK KK KK KK KK KK KK KARA RARA KKK KK KKKKKKKKKKKKKKKK _Start_oversample slave okckckckockckck c kockokck kockckck kokokckck kock REX RAR RARA RAE RARA kockckck kc ck kc ckck kck kck k ck kokock kck kok ck k ko Data output file Outdata slave ovs raw Data output format complex short First line w r t original image 101 Last_line w r t original_image LSS First_pixel w r t original_image 991 Last pixel w r t original image d S Multilookfactor azimuth direction i Multilookfactor_range_direction eS Number of lines multilooked SS Number of pixels multilooked dem First line w r t ovs image Al Last_line w r t ovs_image 133 First pixel w r t ovs image 3961 Last pixel w r t ovs image ADO DER eo Fe Te He ee ON See EERE eoe tee ie Fee See e ae tee ee eee RAR RARA AR AAA te Oo ete ae x End oversample NORMAL kk kk Sk kk Sk kk Sk Sk Sk Sk S kk kk Sk kk Sk kk Sk kk kk Sk kk kk Sk Sk kk kk kk kk Ck Ck Ck Ck Ck A A A A A ko ko KAR 8 3 Algorithm Based on description by Raffaele Nutricato who provided this code In the code look for A RaffaeleNutricato START MODIFICATION SECTION 1 As I explained in the previous e mail range oversampling is obtained as convolution of the zero filled signal with a truncated sinc The loading of the image is performed line by line and co
96. e real part of the file to the phase and the imaginary part to 1 arbitrary 29 3 2 spectral This method is implemented the same as goldstein method w r t the overlap blocksize etc Algorithm per block SIZE SIZE is to perform a 2D FFT of the block and then to multiply pointwise with the kernel which is padded with zeros The kernel is centered around zero frequency 29 3 3 goldstein The algorithm is implemented as 97 e Read in buffer B of PF BLOCKSIZE lines overlap e Get block B B as input block see Fig 29 9 e B fft2d B obtain complex spectrum e A abs B magnitude of spectrum e S smooth A convolution with kernel e S S max S S between 0 and 1 e B B S weight complex spectrum e B ifft2d B result in space domain e f all blocks of buffer done write to disk pixel range blocks _ A A A E i 1 PF_BLOCKSIZE i B i 2 t B1 i E i CEN cr RR E if 2 PF_OVERLAP i O a M PUE EN 8 y D B2 La A ERU T S _ pal complex interferogram ge 7 3 Pd 41 4 1 Per buffer d 4 pee 1 d T 1 bd 1 A 1 x 1 1 Bl 2 2 BI2 nae i 4 1 ox g I T Output for this buffer in this direction i 1 x lx 1 x E da AAA ELLE LLL ee ee ee mum m m m t 4 BI X t 2 PF OVERLAP J 4 N 26t Per block AC nuin nin ania pU E 28 1 1 1 2PF SMOOTH 1 OH 1 i 30 E y 31 Fo I Freee l B11 1
97. e real valued and the computed reference phase is saved to a file Second the reference phase is interpolated to the integer grid of master coordinates A linear interpolation 85 based on a Delaunay triangulation is used The software package Triangle for the Delaunay triangulation is kindly made available by Jonathan Shewchuk Shewchuk 1996 Shewchuk 2002 see also http www cs cmu edu quake triangle html 86 Chapter 27 SUBTRREFDEM In this chapter the processing of step SUBTRREFDEM is described This step requires the steps INTERFERO and COMPREFDEM for obvious reasons In this step the reference phase of a digital elevation model is subtracted from the complex interferogram This is done by complex multiplication with the conjugated as explained in step subtrrefpha chapter 25 27 1 Input Cards SRD_OUT_CINT cint minrefdem raw Filename of output complex interferogram SRD_OFFSET 00 Offset to be applied in azimuth line range pixel direction The synthetic phase image outputed by COMPREFDEM is shifted by the specified offset before subtraction A pos itive shift indicates a shift of the synthetic phase image to the right range up azimuth Example input section E E step subtrrefdem SRD VA Gaunt Outdata cint minrefdem raw SRD_OFFSET qom 27 2 Output Description At successful exit the process control flag is switched on subtr refdem il The output in the products result file looks like
98. e H from the position of M H can be computed by reference needed as shown And then to find theta f B phi r And then find h with the first equation This method is implemented to check our exact method the results are very different In this method the point P does not have to be computed though in order to compute baseline components we will compute a point P on height h evaluate iteratively once for every line A better way might be to use the co registration model to compute the point S 32 3 3 Method schwabisch This method is described in Schw bisch 1995 It is a fast method that yields the radar coded heights It is based on the idea to first compute the reference fase at a number of heights and then to compare the actual phase from the interferogram with these values to determine the height A problem is that the interferogram does not contain the reference phase anymore so that has to be added to the estimated phi It uses a number of steps which are described below 1 Compute reference phase at a NL locations line pixel at NH 23 heights 0 2000 and 4000m for h 0 2000 4000 e ellips a wgs84 a h ellips b wgs84 b h e compute position of master satellite corresponding point P on ellips and position of slave satellite see annex D e Bi ri r2 e r B e store the values and locations and heights Note that the reference phase defined like this typically is something like 5000 rad even for
99. e included the code for a Cooley Turkey algorithm but it handles the data quite slowly compared to an optimized library E 1 Matrix class functions The functions in the matrix class are obtained by bin ctags c types f x matrixbk cut c1 12 50 600 155 0 400 Cc O o0 0c1Qg1C1Qg1C1Q101010 1C i amp BBHKRHHBRHKHPHWAWHWWWAWWAWHWAWNNNNNNNNNNHHHBZBHeeteostetes O coc 100iARcC0rn ocIuoo0 100dc0rnm occiuo0 1o00 mRcnm oc do0 0o0o0 7 cm5 ocKdo0 o0230 crmc oco allocate checkindex clean conj conj correlate diagxmat dotdiv dotmult dumpasc fftshift fliplr flipud getcolumn getdata getrow ifftshift initialize isvector lines matTxmat matrix matrix matrix matrix matxmatT max max mean min min multilook mypow myswap operator operator operator operator operator operator x operator x operator operator operator operator operator operator operator operator operator operator operator operator operator operator lt lt operator operator operator operator gt gt operator operator operator void matrix lt Type gt void matrix lt Type gt allocate uint numlines uint numpixels allocator checkindex uint line uint pixel const void matrix lt Type gt clean sets 2 zero void matrix Type conj matrix lt Type gt conj const matrix lt Type gt 8A matrix lt real4 gt corr
100. e related chapter for details on the magspace and magfft methods 27 Chapter 8 M_OVS In this chapter the processing of step M_OVS is described This step can be run optionally to oversample the cropped data and is done after M CROP and optinally after MLSIMAMP and M TIMING Range oversampling has been implemented by Raffaele Nutricato who uses an oversampling factor of 4 in range for advanced processing in the multi temporal analysis For PS type processing factor two in both directions seems reasonable This avoid aliasing in the spectrum of the interferogram which implies you can interpolate correctly in the interferogram 8 1 Input Cards M OVS OUT master ovs raw Filename of the oversampled data M OVS OUT FORMAT ci2 Output file format M OVS FACT RNG 1 Oversampling factor of output image in range pixels M_OVS_FACT_AZI 1 Oversampling factor of output image in azimuth lines M_OVS_KERNELSIZE 16 Kernel size sinc used for oversampling Example input cards for this step e e comment __OVS___ e M_OVS_OUT Outdata master_ovs raw output filename M_OVS_OUT_FORMAT ci2 Oe U EE ramalts Hy sane ciz cT M_OVS_FACT_RNG 2 range oversampling ratio M_OVS_FACT_AZI 2 azimuth oversampling ratio M OVS KERNELSIZE 16 interpolation kernel length 28 8 2 Output Description The process control flag at the start of the result file is switched to 1 at successful exit oversampl
101. e run a certain step without editing the result file The process control flag in the header has to be reset and the total section in the tail from Start_step to END NORMAL has to be deleted or commented out If the section is not deleted Doris will likely exit but if not the further processing may be affected because wrong values may be used i e read from the result file It is of course possible to change the results parameter values for example correlation value for an estimated offset in the result files so that the altered value is used in the further processing However if you change the strings describing the output Doris will likely protest i e hang or exit 1 3 Outline of this document In Chapter 2 the general purpose cards are described All further chapters describe a certain processing step The Chapter name is equal to the argument of the PROCESS card that should be given in the input file to switch on the processing of that step The step is introduced at the beginning of the chapter In the first section the input cards are described and also an example input file is given In the second section the output is described as well with an example In most chapters there is also a third section describing the implementation Cards that are mandatory are written in boldface sans serif Optional cards are in normal sans serif font Parameters are notated in italic and defaults are underlined Optional parameters are
102. e that for each interferogram of the stack not a large file is created for the master The disadvantage of not filtering the master of course is that a small part of the spectrum of the master is not shared with the slave spectrum yielding coherence loss of the interferogram Further information on the input output of this step can be found in Chapter 15 M FILTAZI 43 Chapter 17 FINE In this chapter the processing of step FINE is described The offset vectors to align the slave image to the master are computed with sub pixel accuracy for a number of locations in the master Over the total image for a large number of windows e g 500 distributed by Doris or from a file with locations in the master coordinate system the offset between master and slave is estimated by computing the correlation of the magnitude images for shifts at pixel level Next in a local neighborhood of the maximum correlation at pixel level these correlations are harmonically oversampled interpolated requires FFT to find the maximum at sub pixel level These offsets are then written to the products result file The offset is computed in the spectral or in the space domain which is implemented to avoid the use of FFT but that is required later anyway and to provide a check of the method in the spectral domain which should be faster The correlation is computed on the magnitude images Though we believe this to be a good method we would like to investigate fir
103. ee GV Pee Wa 13 2 Output Description d e e s ss e seet E m omm ee EO 13 3 Implementation 0 00000000 14 COARSECORR 14 1 Input Cards rmm i e a Ra eee 14 2 Output Description uod sor a a ae e RR we aE ten e ee heen 14 3 Implementation 0 0000000040 14 3 1 Method MagSpace o e 14 3 2 Method magit c sry ss mg 15 M FILTAZI 15 1 Input Caras suit hod xe m Ret dixo po ECC ae DRE Rota 15 2 Output Description n 15 3 Implementation clle 16 S_FILTAZI 17 FINE Tad PUL Cards crate tera estne on ae 2 a me ne tu a CR 6 eet ee ae anc 17 2 Output Description 2 0 o e e 17 3 Implementation e IAS Mags pace s c uu a io ir eE en ad ak aes 17 3 2 OVersample iocos eso ncm uec we wea ce Ride en e Re e Uu 17 3 3 magift 2 220 20 9 en ee at heehee ede ee bee He 18 RELTIMING 18 1 InputCards 2 2 ee ee 18 2 Output Description n 18 3 Implementation celer 19 DEMASSIST TIA Jnput Gards nane xut gs poerem e 9x a NIE mU de da A 19 2 Output Description lt o cs yoy eame vw RPM 19 3 Implementation e 20 COREGPM 20 1 Input Cards i ssa 2s so RR Er mmt 20 2 Output Description lees 20 3 Implementation ees 21 RESAMPLE A A a n a a e E a a e ae A a A a EE E a 21 2 Output Description ss acae si s ai Sioi
104. ee methods Figure 32 4 shows a comparison between schwabisch and ambiguity method It can be seen there is a trend between schwabisch and ambiguity And there is an offset with rodriguez Schwabisch method is always higher then ambiguity suggests error in computation of baseline parameters 114 h m h m 150 100 50 100 150 200 0 200 150 100 50 100 150P 200 250 300 comparison slant to height methods pL q Q P p p ei ambiguity i schwabisch h y thay rodriguez A hod hs I DNE d sl NES M n amp 4 F ae UT 1 1 1 1 1 1 50 100 150 200 250 300 350 400 azimuth line number Figure 32 2 Comparison methods in azimuth direction for line 250 comparison slant to height methods ambiguity schwabisch rodriguez 50 100 150 200 250 300 350 400 450 500 550 600 range pixel number Meal a 2 Figure 32 3 Comparison methods in range direction for pixel 110 115 Table 32 1 Processing with the methods Method cpu options remarks Ambiguity 60 geocoding as well Schwabisch 12 1000 pnts 1d 2 2d 5 Rodriguez z mhap Sy i Sel er Jer bai mi n nd ae Lange Figure 32 4 Comparison ambiguity schwabisch methods for total image Some other tests also showed that the method Schwabisch as implemented in Dori
105. eference phase Only if SRPLDUMPREFPHA SRP_OUT_H2PH filename 80 Request optional output of height to phase factors Example input section for dumping the reference phase E SRPSMETBHOD exact SRP_METHOD polynomial Ge SRELOULIGINGE Outdata cint minrefpha raw SRP_MULTILOOK dd SRP DUMPREFPHA SRP OUT REFPHA Outdata refpha raw step subtrrefpha 25 2 Output Description At successful exit the process control flag is switched on subtrrefpha dl The output in the products result file looks like kk CK kk kk Sk Sk kk Sk Sk Sk KK KK KK KK KK KK KK KK KK KARA RARA RARA kk Sk kk Ck RARA ko ko kokok Start subtrrefpha kk kk kk kk Sk Sk Sk Sk kk Sk S kk Sk kk kk Sk kk Sk kk kk kk Sk Sk Sk Sk Sk kk kk Sk Sk kk Ck Ck Sk Ck Ck A A A A A kc kokok ok Data output file Outdata cint minrefpha raw Data output format complex real4 First line w r t original master 205 Last line w r t original master 14964 First pixel w r t original master 7 Last pixel w r t original master 3998 Multilookfactor azimuth direction 40 Multilookfactor range direction 8 Number of lines multilooked 368 Number of pixels multilooked 499 XX koe ie ee ie ooo Eo o o eoe ok ae ek ak to deat c ec Ecke KS eo ke AA ee ae Xo Kok Xo ko ox End subtr refphase NORMAL Ck Ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck Ck ck ck ck Ck ck ck ck Ck ck ck Ck Ck Ck Ck CK CC CK CC CK CK Sk SS Sk Sk Sk I
106. elate const matrix lt Type gt amp A matrix lt Type gt Mask matrix lt Type gt diagxmat const matrix lt Type gt amp diag const matrix lt Type gt 8B matrix lt Type gt dotdiv const matrix lt Type gt amp A const matrix lt Type gt 8B matrix lt Type gt dotmult const matrix lt Type gt amp A const matrix lt Type gt amp B void dumpasc const char file const matrix lt Type gt amp A friend void fftshift matrix lt Type gt amp A void matrix Type fliplr void matrix lt Type gt flipud matrix lt Type gt matrix Type getcolumn uint pixel const matrix Type matrix Type getdata window win const matrix Type matrix Type getrow uint line const friend void ifftshift matrix lt Type gt amp A void matrix lt Type gt initialize uint numlines uint numpixels bool matrix Type isvector const uint matrix lt Type gt lines const return number of lines matrix lt Type gt matTxmat const matrix lt Type gt amp A const matrix lt Type gt amp B matrix lt Type gt matrix constructor 0 arg matrix Type matrix uint lines uint pixels matrix Type matrix const matrix lt Type gt amp A matrix Type matrix window win const matrix lt Type gt amp A matrix lt Type gt matxmatT const matrix lt Type gt amp A const matrix lt Type gt amp B Type max const matrix lt Type gt 8A Type max const matrix lt Type gt 8A uint amp line uint amp pixel real8 mean
107. elation windows in the conven tional method are therefore prevented Furthermore the coregistration improves in case of large base lines and strong topography The improvement is especially significant in case of X band data An analy sis of the performance of the implemented DEM assisted coregistration for various sensors is described in Arikan et al 2008 Nitti et al 2008 19 1 Input Cards DAC IN DEM filename filename of input DEM gtopo30 File is assumed to be stored in a raster Major row order from North to South line by line See also internet links at Doris home page for available DEMs DAC IN FORMAT 12 12 BIGENDIAN R4 R8 format of input DEM on file signed short for gtopo30 or real4 or real8 input ma trix is raw binary data w o header endianness of host platform is assumed except for I2 BIGENDIAN DAC IN SIZE 6000 4800 Number of rows and columns of input DEM file Default is set to tile w020n90 DEM DAC IN DELTA 0 00833333333333333333 deltalon Grid spacing of input DEM in decimal degrees latitude longitude Default is equal gridspacing default set to tile w020n90 DEM DAC_IN_UL 89 995833333333333 19 995833333333333333333 Coordinates of UL upperleft corner in decimal degrees latitude 90 90 longitude 180 180 Default is set to tile w020n90 DEM It is interpreted as max latitude min longitude in source 50 DAC IN NODATA 9999 Identifier to ignore data in input DEM with this value De
108. equation and the following can be easily verified by substitution of the values for the first last line pixel line number The azimuth time of a certain line number ta is computed as 1 tape D 32 l i F PRF Where ta is the azimuth time to line 1 first line And the line number given a certain azimuth time ta can be computed as 1 1 PRF ta ta1 D 33 Where tai is the azimuth time to line 1 first line 149 Doppler centroid The Doppler centroid frequency azimuth is computed as a second degree polynomial foo ao ea zo 02 D 34 p 2 RSR fear Where pis the pixel number starting at 0 a is read from the leader file In the master result file it are the variables e g Xtrack_f_DC_constant Hz early edge 117 3210000 Xtrack_ f_DC _linear Hz s early edge 72338 0000000 Xtrack_f_DC_quadratic Hz s s early edge 455000000 00000 This frequency is used in the azimuth filtering and in the resampling It should also be used if the complex SLC data is harmonically oversampled as is done in the range filtering routine but we did not implement this yet But for fpc 150 Hz this should not have any effect assuming PRF ABW 300Hz Since for a signal f t with Fourier transform F w FT f t lt gt Fu D 35 dht El Ffw wo D 36 D 37 The azimuth spectrum of a SLC image processed on a certain Doppler frequency fpc spectrum shifted to this frequency can be s
109. er result file after the step readfiles at place number of lines original of datafile A function normalize is called to do the normalization so it is easy to change the implementation to a different normalization It has been noticed that for higher order polynomials the normalization factor is very important to obtain a stable estimate 148 D 4 2 Evaluation of polynomials Evaluation of the polynomials should be done by normalizing the data as indicated above Something like const real8 const real8 minL master originalwindow linelo const real8 maxL master originalwindow linehi const real8 minL master originalwindow pixlo const real8 maxL master originalwindow pixhi matrixreal4 N 1 I axis linenumbers normalize IN axis matrixreal4 f polyval IlN axis p _axis coefficients degree NOTE It is faster to evaluate a polynomial on a grid than point by point D 5 SAR System parameters D 5 1 Azimuth PRF The actual pulse repetition frequency PRF is computed based on the data in the SLC leader file However the actual value as read from the leader file is used after private communications with ESA helpdesk It is defined as NI 1 dta PRF D 31 Where PRF is the pulse repetition frequency in Hz Nlis the total number of lines lastline firstline dt is the azimuth time of the last line minus the azimuth time of the first line or the acquisition time of the image This
110. et to BLOCKSIZE 2 1 the maximum value for this parameter then each output pixel is filtered based on the spectrum that is centered around it This is probably the best but may be time consuming PF BLOCKSIZE 32 This card is for method goldstein and spectral only Size of the blocks that are filtered This must be a power of 2 It should be large enough so that the spectrum can be estimated and small enough that it contains a peak frequency 1 trend in phase 32 is recommended PF KERNEL 3 1 8 1 8 1 3 This card is for method goldstein and spatialconv only 1D Kernel function to perform convolution First the number of elements in the kernel is given then the values The kernel is always normalized to 1 by dividing the kernel by the sum of the absolute values of the kernel For method goldstein defaults to kernel 1 2 3 2 1 This kernel is used to smooth the amplitude of the spectrum of the complex interferogram The spectrum is later scaled by the smoothed spectrum to the power alpha For method spatialconv Default is a 3 point moving average 1 1 1 convolution The real and imaginary part is averaged seperately this way For more info see implementation section The output matrix has a zero valued edge of size floor kernel 2 PF IN KERNEL2D filename This card is for method spatialconv and spectral only Name of ascii input file to specify a 2D spatial kernel function This file must start with a 1 line header containing numlines num
111. f course the main advantage is a factor 2 reduction in the size of the output file Now an example follows for the error in amplitude and phase for a complex value of about 100 100 If the actual interpolated complex value equals 100 5 100 5 then the error in the magnitude approximately is em 100 4 100 100 4 100 5 100 52 100 5 100 52 0 5 21 1 If the actual value did equal 100 5 100 0 then the error in the phase is approximately ep 100 arctan 100 100 arctan 100 5 100 arctan 100 5 100 0 3 21 2 These are worst case scenarios If the complex value is larger then the relative error decreases Note that the maximum for a signed short integer is 21 32768 64 21 3 2 Interpolation Kernels In this section the available kernels are defined See also Hanssen and Bamler 1999 The KNAB interpo lation kernel is described in a IEEE letter of 2003 The Raised Cosine interpolation kernel is described in a article in J Of electromagnetic waves 2005 Cho et al sin Tx sinc r TX 0 z 0 5 rect z lt 0 5 a 2 0 5 1l a lt 0 5 0 lx gt 1 ilz tri z 1 z e lt 1 a 1 a 2 2 a 3 2 1 0x lx lt 1 i z o z 5alz 8a z 4a 1 lt 2 lt 2 0 2 lt lx a 5 B 5 o B 2 z a 8 3 lz 1 TE alx 5a 8 a 8a 38 2 4a a 88 x 218 x 188 L 6 8 16
112. f step COMPREFPHA is described This step can be performed as soon as the precise orbits are known The recommended approach is to compute this only after the computation of the interferogram and then use the step SUBTRREFPHA to subtract it This step is not required if method exact is used in step SUBTRREFPHA The flatearth correction the phase caused by the reference surface WGS84 for now is computed in this step For a certain line pixel in the master image the corresponding coordinates x y z of the master and slave satellite and the point P on the reference ellipsoid are computed utilizing the set of equations Doppler range and ellipsoid equation see annex D Then the parallel baseline and the phase is computed The parallel baseline B is defined as M S are positions of master slave P is the position of the point on the reference surface By d M P d S P 24 1 The phase of a pixel in the master image is defined as 4 Eur o M P 24 2 The reference phase for this pixel is defined as Bi 24 3 The reference phase is computed in a number of points distributed over the total image in this way after which a 2d polynomial is estimated least squares fitting these observations So a plane can be fitted by setting the degree to 1 A 2d polynomial is defined as d i f z y p3 y oj 4 2 Hy 24 4 i 0 j 0 Thus the order of the coefficients line pixel is degree d d 0 Aoo
113. fault is set to tile w020n90 DEM DAC OUT DEM filename Request optional debug output to float file of input DEM per buffer cut to the interfero gram window Info on these files is written as DEBUG DAC OUT DEMI filename Request optional debug output to float file of interpolated input DEM cut to the interfer ogram window Info on these files is written as DEBUG DAC OUT DEM LP demheight lp raw Filename of output DEM height in radar coordinates Example input section e comment DEMASSIST____ c DAC IN DEM final wanaF2835 dem DAC IN FORMAT r4 DAC_IN_SIZE 3601 3601 DAC_IN_DELTA 0 000833333 0 000833333 DAC_IN_UL 40 27 the center cn of UL corner p DAC_IN_NODATA 32768 DAC_OUT_DEM dem_dac raw c DAC_OUT_DEMI demi_dac raw c DAC_OUT_DEM_LP demLP_dac raw 19 2 Output Description At successful exit the process control flag is switched on demassist 1 The output in the products result file looks like KKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKK Start dem assist c ck c ck ck ck ck 0k ck ck KK ck ck ck ck ck Ck ck Ck Ck ck c ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck koc KKK ck ck ko KKK ck Sk ck ko ck ko kx ko ko kx x DEM source file final wanaF2835 dem Min of input DEM 92 Max of input DEM 1815 First line w r t original master 3053 Last line w r t original master 8052 First pixel w r t original master 1714 Last pixel w r t origin
114. fo pair see chapter 2 Filename of the master of the topo pair To obtain the orbit and other parameters for the topo master DLIN TOPOSLAVE result file name Filename of the slave of the topo pair To obtain the orbit and other parameters for the primary slave DI IN TOPOINT result file name 102 Filename of the interferogram result file of the topo pair To obtain the name and di mensions of the unwrapped topography interferogram DI OUT SCALED filename Filename for optional debug output of a real4 file with scaled with ratio of perpendic ular baselines unwrapped topo interferogram Example of the cards for this step c c comment DNS ASAS SMN E DI_OUT_FILE Outdata difg raw c DI IN TOPOMASTER data project topo master res if 4 pass method DI IN TOPOSLAVE data project topo slave res DI IN TOPOINT data project topo products res CED TECUTES CATED Outdata scaled raw debug 31 2 Output Description In the defo processing result file for the products the process control flag for dinsar is switched on dinsar JL A complex real4 mph file is created with the wrapped differential phase The amplitude is the same as that of the original deformation interferogram A complex value 0 0 indicates unwrapping was not ok If the the debug version of Doris is used compiled with DEBUG then ascii matrices are dumped for Linenumber Pixelnumber Bperptopo Bperpdefo and Ratio Fig
115. g this line This may have to do with my implementation of complext h from GNU Experience Running Doris Mohanty KK 132 B 5 4 Some notes on installation on Linux X86 Since July 2000 Doris can be installed to Linux X86 systems The code needed some changes Please refer to the general problems section as well e Byte order on X86 systems Use the functions htonl and htons to read the record length and data of the SLC volume leader and data file e File io gcc compiler only ios ate does not act as expected It should open an input file stream at the end but it does not ios app ofstream class does work as does a file seekg 0 ios end ios in or out has to be set even for statements like ifstream ifile ios nocreate The default and logical that a file is an input stream if it is declared ifstream is not true for gcc compiler e strcmp The statement stremp word V0 gives a memory fault if tested B 5 5 Some notes on installation on Window running Cygwin Bert Kampes reported installing Doris on Windows NT and XP running Cygwin without trouble beginning of 2002 The run installation scripts needed some small changes Here is the summary It is assumed that a full version of Cygwin installed including tcsh developing tools e There is no ksh and csh for Cygwin in standard setup Therefor change the first line of the configure script to tcsh and use the new versions of the helpdori
116. h 0 Therefor the reference phase for h 0 is set to 0 because in the unwrapped interferogram the reference phase is removed if the phase of the unwrapped interferogram is 0 then this should yield a height of 0 and the reference phase for height h is set to refphaseh refphase0 This makes the computations as done later a little stupid to estimate coefficients which are by definition equal to 0 An other possibility is not to do this here but later when the functions are evaluated to add the reference phase to each pixel have tested this and the results are identical 4m 2 Compute for each location a polynomial 1d degree 1dD NH 1 to describe the height as a function of reference phase at these points For each location it holds 6 for height i h 0 00 01 a 32 24 h 2000 ag ai 029 32 25 h 4000 ag ai 2 0242 32 26 So it is easy to solve exact for a per location 113 3 Compute 1dD 1 NH polynomials to describe the coefficients of the previous step as a function of location now for a random location the coefficients of height as a function of the reference phase can be computed For example for ag computed at NL locations l p a 2d polynomial can be used Og Soo Biol borp 8200 Soe ul 32 27 A linear system can be easily set up and solved least squares by cholesky factorization A rescaling needs to be applied to avoid instability The system can be solved simult
117. hase One can use an awk like to make a grid awk BEGIN for i 100 i lt 25200 i i 500 N for 3 750 3 lt 5400 3 3 200 igwexbewgxE Sa Bal Vin pal a exui A OS MO MS and a card in the input file BEESSENSPOS positions cz abit This step used to be named FlatEarth This explains the prepended FE_ s instead of some abbreviation for reference phase 78 Example input section s comment COMPREFPHA E FE_METHOD porbits FE_DEGREE E FE_NPOINTS 201 24 2 Output Description At successful exit the process control flag is switched on comp_refpha 1 The output in the products result file looks like Wood A oe ie ae ee ee Fe eK kt AAA Fe He He ae tele ee AIR ae eoe AR A ie ke te toe ee es Hee Se he KIRKER He RER Xo es ee oe Start comprefpha obckck coke oe ck eoe ko ck e ok cec ke ke kc kc ACER RARA AREA AAA AAA ck kok kock kokck kc k kockok k kk kok k kk Degree flat 5 Estimated coefficients flatearth 5 17144173e 03 Qu O 4 03705656e 03 2 17736976e 01 2 05452064e 06 2 15880157e 07 1 27934869e 05 8 80499980e 10 i 643303 e TN 8 28354289e 11 6 04376860e 10 1 64239900e 13 A ADOS dL 1 48243991e 14 5 52622526e 14 IAS SON SITE il 2 11724810e 17 74597475e 20 2 64743817e 18 4 10973605e 18 3 09581184e 17 6 94891272e 16 Dk ck Ck ck ck ck ck ck ck Ck Ck ck kc ck ck ck Ck Ck ck ck ck ck ck kk ck ck ck ck ck AAA ck ck ck ck koc ck ck ck kk ck c
118. hase is actually computed as I M S cos sin p D 20 because the complex conjugated of the reference phase equals ae a cos sin cos sing D 21 The complex coherence between two images is defined as see Touzi et al 1996 B E M S YE EMMA FiS SS Where Ef is the expectation is the complex conjugated Ye is the complex coherence Mis the complex master image S is the complex slave image possibly minus complex reference phase S S R The coherence is defined by and its estimator as D 22 D x pum 5 9 The correlation between two images is defined by see B hr and V gtle 1991 p M S _ E M S E M E s 0 24 var M var S J E M M E M E M E S S E SHE S Thus the mean is first subtracted in comparison to the coherence The coherence is equal to the correlation only if E M E S 0 D 23 MM 2 A problem is that the estimator for the coherence and correlation is biased For small window sizes its outcome is too high This probably also causes the problems in the coarse coregistration where the most likely offset is not selected based on its correlation value but on its consistency 147 D 4 Polynomials A 1d polynomial is defined as f t Y aya D 25 A 2d polynomial is defined as d i fw 5 5 o jz Jy D 26 i 0 j 0 Thus the order of the coefficients line pixel is in
119. he routine splineinterpol file utilities c where the coefficients are computed NATURALSPLINE is de fined This sets the boundary condition to use zero second derivative at the borders Otherwise the first derivative is set to a specified value This does not seem to make a big difference 19 Chapter 5 M_CROP In this chapter the processing of step M_CROP is described This step normally is the third one that is run after M READFILES and M_ORBITS It requires the SLC data file on disk or cdrom For ENVISAT a utility is called that does the work In this step the SLC datafile is put on disk in a raw pixel interleaved 2x2byte signed short integer complex format The reason for this step is that we normally work with the SLC images on cdrom and that we want to have the files on disk to perform operation requiring both the images It also serves as a common format for different input A few checks are performed regarding the number of lines which is written in the header of the SLC data file as well as in the leader file The image is read written line by line no data conversion takes place though a cutout can be made If you are working on a little endian platform X86 PC then the data is converted from big endian which is the CEOS format 5 1 Input Cards M_CROP_IN filename Filename of the SLC data file M CROP OUT master raw Filename of the raw data output file M DBOW linelow linehi pixellow pixelhi Master ou
120. hen input is not copied if no parameter is given it defaults to ON do copy MEMORY 500 With this card the user can indicate the maximum amount of memory to be used by the processor in Megabytes It is advised setting this lower than the maximum available amount because it may be somewhat inaccurate up to a factor 2 particularly due to temporary copies created by the copy constructor A lot of routines actually try to use a minimum of memory even if this card is set to a large value PROCESS M READFILES M PORBITS M CROP M SIMAMP M TIMING M OVS M FILTAZI FILTRANGE S READFILES S PORBITS S CROP S OVS S FILTAZI COARSEORB COARSECORR FINE RELTIMING DEMASSIST COREGPM RESAMPLE INTERFERO COMPREFPHA SUBTRREFPHA COMPREFDEM SUBTRREFDEM COHERENCE FILTPHASE DINSAR UNWRAP SLANT2H GEOCODE With this card the processing steps that have to be processed can be switched on More than one PROCESS card can be specified in the input file An ONLYPROCESS card overrides possible PROCESS cards Description of these steps can be found in the introduction chapter 1 and in the following chapters At least one PROCESS or ONLYPROCESS card is mandatory ONLYPROCESS same arguments as PROCESS card With this card a processing step that has to be processed can be switched on Overrides possible PROCESS cards This card also automatically switches BATCH ON At least one PROCESS or an ONLYPROCESS card is manda
121. her with the height width in pixels This card is especially useful in combination with the M DBOW GEO and S_DBOW_GEO cards see Chapters 5 and 11 RS_SHIFTAZI ON OFF If ON then it is accounted for non centered azimuth spectrum of the data The azimuth interpolation kernel is shifted to the Doppler center frequency before resampling If the fDC is small users could switch this card to OFF Example input cards for this step E E comment ___RESAMPLING SLAVE E c RS METHOD cc4p RS_METHOD cc p c RS METHOD ts6p c RS_METHOD ts8p c RS METHOD ic M Io RSHOU IEA Output 01393 resampled RS_OUT_FORMAT 22 c RS_DBOW IGGL 2105 SOL 700 RS_DBOW_GEO 59251125 So VSZ SOOO 1101010 21 2 Output Description The process control flag at the start of the slave result file is switched to 1 at successful exit resample il Example of output of this step in slave result file KC ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck Ck ck ck ck Ck ck ck Ck Ck Ck ck Ck Ck Ck Ck C C KKK KK CC CK CIS Sk ke Sk Sk I Sk ke ke ke kk ko kc kc kc KKK x Start resample kckCckckckckckckckckckckckckckckckck ck ckckckck ck ckckckckckckckckckck ck KK KK ckckckckckckckckckckckckck ck ck ck ck ck ck ck ckckckck k ck kk Data output file Output 01393 resampled Data output format complex short Interpolation kernal 6 point cubic convolution First line w r t original master LOCI Last line w r t original
122. hifted back to zero by multiplication in the space domain by the term e 327 LS line D 38 See also any signals and systems book or e g Geudtner 1996 The spectrum can be shifted back to the original doppler centroid frequency by multiplication by after e g interpolation foc Pag line D 39 e 7 Proper care should be taken to get the correct line number in both situations D 5 2 Range RSR The range sampling rate RSR is defined as Np 1 r RSR 0 001 D 40 Where RSRis the range sampling rate in MHz Np is the number of range pixels dt is the zero Doppler two way time to the last pixel minus the range time to the first pixel in milliseconds 150 pixel number The pixel number p 1 Np given a certain one way range time ta can be computed as p 14RSR 2 t tri D 41 Where t is also one way The one way range time for a given pixel number p can be computed as p 1 os 1 RSR D 42 Where RSR is in Hz tai is the range time to pixel 1 first pixel in seconds The range is of course equal to PSA ee D 43 Where cis the speed of light constants h 299792458 m s D 6 Doppler range and ellipsoid equations The following three equations are used regularly throughout the software to compute the point P that corre sponds to a certain line and pixel in the master or slave image see also Geudtner 1996 Precise orbits are necessary 1 Dop
123. his card is not required COMMENT After this card everything up to a newline is ignored A space after this card is not required SCREEN DEBUG INFO PROGRESS WARNING ERROR This card controls the level of standard output It is recommended to start with this card since it is in effect only after it is read BEEP OFF WARNING ERROR PROGRESS ON This card controls the level of beeping BATCH ON OFF Specifies to run the processor in non interactive mode If this card is omitted then the processing is done in interactive mode asking to press a key before each step BATCH can be specified without arguments which means non interactive processing If there is an ONLYPROCESS card present this forces BATCH ON 10 OVERWRITE ON OFF Specifies whether or not to overwrite existing files If this card is omitted files are not overwritten if no parameter is given it defaults to ON do overwrite PREVIEW OFF ON XV Specifies whether or not to generate SUNraster preview files with the help of the utility program cpxfiddle download and install seperately from Doris website Default is OFF since this program may not be installed If ON shell scripts are created in the working directroy which create the SUNraster file if run If XV is given also the command is given to view the generated file with xv LISTINPUT ON OFF Specifies if the input file has to be copied to the logfile If this card is omitted t
124. hive there are a number of utilities included which are required for optimal processing They can be called from within Doris The utilities in the bin directory are explained in section B 6 If these files are in your path and they can be executed then they are installed ok Doris writes the calling syntax of these utilities to standard out as INFO so you can repeat the commands For example assuming the run script is used grep plotoffsets Outinfo out grep plotcpm Outinfo out 126 yields for example Outinfo out input fine cpm 4926 INFO plotcpm CPM Data 1 5000 1 1000 amp This command can be repeated from the prompt B 2 Additional programs To obtain a full version of Doris the getorb precise orbits GMT mapping tool helper for fine co registration and gv or ghostview to display postscript files in csh script plotscript should be installed on your system These programs can be obtained freely but are not included in the Doris distribution At our homepage you can find out more on how to obtain these packages http enterprise r tudelft nl doris B 3 Running the Doris software You can run the Doris software by making an input file as described in this user s manual For full function ality make sure that getorb ghostview gmt are in your path add them in your home cshrc file or in the Shome login file The command line options for Doris are doris v Return version number doris h search pat
125. how this is done If the VECLIB library is not used slower internal functions for FFT and matrix multiplication are used If you like to use your own FFT only the function four1d has to be changed since the 2d function call this one sequentially If the LAPACK library is not used slower internal functions for cholesky are used The matrix class is used as a container class for part of the images It has not been defined as a class radarimage because in that case it would have been difficult to perform operations on the images if they didn t fit in the memory as a whole However it might be a good idea to define functions such as phasefilter for the matrix class This would result in calls like matrix lt complex real4 gt BUFFERMASTER Container for i 0 i lt NUMBUFFERS i Read in phase image BUFFERMASTER readfromfile filename windowsize formatflag filter this buffer BUFFERMASTER phasefilter parameters Write to outputfile BUFFERMASTER writetofile filename2 formatflag which seems very readable and maintainable only the member function phasefilter of the matrix class has to be changed if something has to be added If you do not have access to LAPACK but you have a different library we advice you to use that one The Cholesky factorization as implemented internally has not been optimized in any way The same holds for the VECLIB library particularly the FFT routines We hav
126. ic aperture radar systems and signal processing John Wiley amp Sons Inc New York Eineder 2003 Eineder M 2003 Efficient simulation of SAR interferograms of large areas and of rugged terrain IEEE Transactions on Geoscience and Remote Sensing 41 6 1415 1427 Gatelli et al 1994 Gatelli F Monti Guarnieri A Parizzi F Pasquali P Prati C and Rocca F 1994 The wavenumber shift in SAR interferometry EEE Transactions on Geoscience and Remote Sensing 32 4 855 865 Geudtner 1996 Geudtner D 1996 The interferometric processing of ERS 1 SAR data Technical Report ESA TT 1341 European Space Agency Translation of DLR FB 95 28 Geudtner and Schwabisch 1996 Geudtner D and Schwabisch M 1996 An algorithm for precise recon struction of InSAR imaging geometry Application to flat earth phase removal phase to height conversion and geocoding of InSAR derived DEMs In European Conference on Synthetic Aperture Radar K nigswin ter Germany 26 28 March 1996 K nogswinter Germany Ghiglia and Pritt 1998 Ghiglia D C and Pritt M D 1998 Two dimensional phase unwrapping theory algorithms and software John Wiley amp Sons Inc New York Goldstein and Werner 1998 Goldstein R M and Werner C L 1998 Radar interferogram filtering for geophysical applications Geophysical Research Letters 25 21 4035 4038 Hanssen and Bamler 1999 Hanssen R and Bamler R 1999 Evaluation
127. ike to stay independent from a particular orbit format and these computations can be done fast anyhow The velocity can be interpolated by the derivative of the piecewise polynomials see source code or numerical recipes The accerelation can be interpolated by the second derivative of the piecewise polynomials see source code or numerical recipes If less points are known than the typically 21 of getorb one wants to use the SLC datapoints for a quick look analysis for example then this kind of interpolation probably does not work very well In future we will include an option to interpolate by a low degree polynomial which is estimated least squares from the datapoints Getting the derivates at any point is straightforward in this case As a satellite moves very smoothly a polynomial of a lower degree might even be nearer to the true orbit then a piecewice polynomial In future we want to model the baseline Bh Bv as a function of azimuth time by a first order polynomial This probably is more efficient than computing the positions of the sensors each time the baseline is required We do not know what the best way is to do this D 8 Format of the products Start at azimuth line 1 range pixel 1 near range Data is written line by line major row order We give the binary data a raw extension The complex interferogram is written pixel interleaved see D 2 Each complex pixel is written as 4 byte real 4B imagina
128. in between square brackets The definitions used in the software are described in Annex D Here amongst others the baseline definition the file formats used and normalization of polynomials are described In Annex B the instalation of Doris is described It also contains a small trouble shoot section Annex C shortly describe third party packages that should be installed for a complete version of the Doris software Also some utilities we have developed are described The matrix class which comes with the Doris software is described in Annex E This matrix class can be freely used in other non commercial programs Finally Annex F describes how to add a module to the Doris software Extention of the software is encouraged Chapter 2 General Cards This chapter deals with input cards that are not specific for a certain processing step the general input cards For example such a card could specify whether you like to do batch processing or interactive processing These cards are best placed at the start of the input file They do not generate any specific output An example of the header of the input file is given in section 2 2 2 1 General Input Cards After this card everything up to a newline is ignored A space after is not required i leal After this card everything up to a newline is ignored A space after is not required C After this card everything up to a newline is ignored A space after t
129. in one file 153 Table D 3 Format of a hgt height file unwrapped complex interferogram and others Actually is a major row order band interleaved data with 2 float 4B canals amplitude phase 1st pixel 2nd pixel 1st pixel 2nd pixel Tstline Amplitude Amplitude Phase Phase 2nd ine Amplitude Amplitude Phase Phase 3rd line Amplitude Amplitude gt gt EET NN ERE ORE I ME RS ER inel Amplitude Amplitude Phase Phase Computations are done in general in the original master system no matter if cut out or multilooked Time system orbit is in seconds of day Ephemerides orbit system is more or less WGS84 Matrix class A matrix is starts at 0 0 etc Offset Offset for a certain point is defined as coordinate in slave system coordinate in mastr system offsets 154 Annex E Matrix class The template matrix class called matrixbk that is provided with the Doris software can be used for other applications as well Please refer to the file matrixbk_test c for an example how to use this class in your own programs do not claim it is the best fastest implementation ever i just find it very useful and the routines are checked and working fine The Makefile shows how to compile it Also see Annex B The data has to be linear in memory for the VECLIB library thus this has been done and it used in some other functions to speed them up see the constructor
130. include refdem Me OEC DURS ENG S To Te e I RIS SI complex image CIOJSL AUT CYOJBI Output coh raw real COH MULTILOOK 10 2 COH WINSIZE 10 2 90 28 2 Output Description At successful exit the process control flag is switched on coherence il Example output section for this step in products result file ck c 0k ck ck Ck ck kk ck Ck ck KK KK KKK KKK ARA ck ck ck ck ck ck ck RARA ck ck ck AAA RARA ck kk ck AAA ko ko ck kv kv ko ok ok x Start coherence Ck ck ck ck ck ck Ck Sk ck kk KK KKK KKK KKK KK KKK KKK KKK ck Ck ck KKK kk ko kk kk ko kk kk ko kk Sk kc ko KKK KKK KKK ethod INCLUDE REFDEM Data output file AAO meen 35 CON Data output format real4 First line w r t original master 3060 Last line w r t original master 8052 First pixel w r t original master 1719 Last pixel w r t original master 2 LO Multilookfactor azimuth direction 5 Multilookfactor_range_direction T Number of lines multilooked 998 Number of pixels multilooked SO ck ck 0k ck ck Ck ck kk kk ck kk ck ck ck ck ck ck ck ck ck kc KKK ck ck ck ck ck ck ck ck ck ck ck ck ck ck RARA ck ck ock ck Sk ko ko ko ck kv kv ko ok ok End coherence NORMAL Ck ck ck ck ck ck KKK KKK ck ck ck Ck kk kk ck kk KKK KKK KKK KKK ck ck AAA RARA ko kk kk ck kk Ck kc ko kk ko ko kk ck ko The output data file must be viewed with an external package like Matlab for now 28 3 Implementation The images are read in buffers for
131. ined system of 3 equations yields the next solution dx and the new values for the unknowns become dz dzo dx which are used to compute dy and d A The solution is updated until convergence Az 1e 6 meters Where dy contains the observations set of equations dx contains the unknowns coordinates of P dA contains the partials evaluated for previous solution To solve for the azimuth time if the coordinates of a point on the ground is known only the Doppler D 45 equation needs to be used the derivative with respect to azimuth time of this equation equals OB gap ote dx ii D 51 The solution is equal to use approximate solution tag to evaluate these expressions E ta D 52 Ola and taz tao taz D 53 The solution is updated until convergence At lt 1e 10 seconds The range time is then computed as in equation D 46 fae MESA D 54 C D 7 Orbit interpolation We assume the precise orbits are given some time before the first and after the last azimuth line Normally we use getorb to obtain satellite ephemerides with a time interval of 1 second approximately 21 datapoints for a frame of typical 15 seconds Natural cubic splines are then used to interpolate the orbit Because these splines do not behave very well at the edges we use some points before after the first last line Note that the x y and z coordinate are interpolated independently The Delft precise orbits and the getor
132. ing their correlation to 0 00001 in the output section of the FINE processing step Also the observations itself and some statistics are plotted w tests a large value indicates an unreliable estimate The script can be adapted to your own wishes it simply calls GMT see Wessel and Smith 1998 based on the ascii data file CPM_DATA This step is important since the interferogram is sensestive to mis alignments of slave on master There fore we always took a very cautious approach However that meant running this step editing the result file running again etc which got quite cumbersome To reduce the manual effort we introduced a card CPM_MAXITER that performs a number of iterations automatically It also should remove no more windows than necessary for a good fit have experimented with values like 20 for this card having say 600 windows after step FINE If you want to approach that after each computation you want to have full control what to do simply set this card to 0 The first run of coregpm for the area of Fig 23 1 is shown in Figures 20 1 20 2 and 20 3 We have used a polynomial of degree 1 and a threshold of 0 4 here The second run of coregpm is shown in Figures 20 4 20 5 and 20 6 In the products result file the outliers are artificially set to O correlation thus being below the threshold to exclude them from the least squares estimation After this run we continued with the resampling Degree d 1 is enough
133. inline uint filesize 87 getoverlap ioroutines c window getoverlap 88 getoverlap ioroutines c window getoverlap 89 readcoeff ioroutines c matrix lt real8 gt readcoeff 90 fillproductinfo ioroutines c void fillproductinfo 91 assert ioroutines c void assert 92 assert ioroutines c void assert 93 tolower ioroutines c void tolower char xs 94 toupper ioroutines c void toupper char xs 95 printWARNING ioroutines h inline void printWARNING 96 ERROR ioroutines h inline void ERROR char ch ONE27 97 ERROR ioroutines h inline void ERROR const charx file int32 line char ch ONE27 98 WARNING ioroutines h inline void WARNING char ch ONE27 99 PROGRESS ioroutines h inline void PROGRESS char ch ONE27 100 INFO ioroutines h inline void INFO char ch ONE27 101 DEBUG ioroutines h inline void DEBUG char ch ONE27 102 DEBUG ioroutines h inline void DEBUG const char file int32 line char ch ONE27 103 initialize orbitbk cc void orbit initialize const char file 104 computecoefficients orbitbk cc void orbit computecoefficients 105 getklokhi orbitbk cc void orbit getklokhi real8 t 106 getxyz orbitbk cc cn orbit getxyz 107 getxyzdot orbitbk cc cn orbit getxyzdot 108 getxyzddot orbitbk cc cn orbit getxyzddot 109 Ip2xyz orbitbk cc int32 Ip2xyz 110 xyz2orb orbitbk cc int82 xyz2orb 111 xyz2t orbitbk cc int32 xyz2t 112 xyz2lp orbitbk cc int32 xyz2lp 113 ell2lp orbitbk cc int32 ell2lp 114 Ip
134. is to use first lines second pixels e g for the order of input arguments The line direction azimuth corresponds to the vertical y The pixel direction range corresponds to the hori zontal x Note that other software may use x before y e The first line pixel of an image is indexed as 1 this may be a bit unusual In the software the first index of an image in a matrix is equal to 0 To index a matrix use MAT y x i e as in linear algebra Matlab etc e The name format and dimensions of the current master slave interferogram are stored in information structs that are filled by reading the corresponding result files The files do not have a header e Generally all coordinates are in the master radar coordinate system The first and last line are given as well as multilook factors for both directions If for example an interferogram is multilooked at the generation with a certain factor and later on again at the subtraction of the reference phase then the first line is still the same see Figure 1 3 e In our view the output files can only contain the results of one algorithm per step So it is not possible to re run a step for example with another algorithm without deleting the previous result New users may be confused by this approach but Doris should give an appropriate error message if it is attempted to run a step twice e Temporary files are created during the processing Their names always start with scratch
135. k 2002 see also http www cs cmu edu quake triangle html The finally optained master slave offsets per master pixel are saved to a file and used in the RESAMPLE step to resample the slave image to the master geometry 52 Chapter 20 COREGPM This chapter describes the processing step COREGPM is described coregistration parameters computation of a polynomial that models the alignment of slave on master Based on the estimated offsets computed in step FINE a 2d polynomial model of certain degree of the coregistration is computed When the step DEMASSIST is applied the polynomial is estimated through the residuals between the FINE and DEMASSIST results This appears necessary because we experienced a remaining trend in the residuals In this case a 1 degree polynomial seems sufficient A least squares solution is used solution by cholesky decomposition of the normal matrix Data may be excluded a priori by setting a threshold value for the correlation Data can also be excluded by editing the products result file and artificially decrease the correlation for a certain offset window After the computations the residuals between the estimated model and the observed offsets are plotted with the csh script plotcpm These plots are useful to iteratively come to a good transformation model changing CPM THRESHOLD or CPM DEGREE for each iteration or identify and remove some estimated offsets the observations blunders by sett
136. k Ck Ck A A A Mk A kc kokokok The time is in seconds of day which can be computed as fractional day 60 60 24 or hours 60 60 min 60 Sec If card M ORB DUMP then an ascii file masterorbit dat is written with the computed t x y z xdot ydot zdot xddot yddot zddot 4 3 Implementation Based on the UTC time of the image and the PRF and number of lines basically the program getorb is called through a UNIX system call Functions from the standard library ctime or time h are used for the time conversions This call is echoed to the screen as DEBUG The commands can be executed stand alone as well oe getodr 950727094923 950727094945 1 data orbits ERS2 ARCS dummyout oe oe read in ODR file name from dummyout ODR 422 untar data orbits ERS1 ARCS ODR 422 getorb 950727094923 950727094945 1 data delftorbits ERS1 gt siciealicClooucloyiic The ephemerides are first written to a dummy file named scratchorbit and later placed in the result file without the velocity output only t x y z It has been noted that sometimes this scratch file is not automatically removed This file can be savely removed by hand Natural cubic splines are used for the interpolation so it may be wise to have a short time interval and some data before the first and after the last line If there is a section with the ephemerides of the SLC leaderfile in the master result file then this section is removed In t
137. k ck ck ck ck ck ck ck ck Ck Ck Ck Ck Ck Ck CC C C C C C CC CC CIC CCCII Sk Sk Sk S ke Sk Sk Sk Sk kc k kc kc kc AR Data output file Output cint raw Data output format complex real4 First line w r t original master VO Last_line w r t original_master eu AA Gas a OOS abso E E 501 Last_pixel w r t original_master 700 Multilookfactor_azimuth_direction RO Multilookfactor_range_direction 2 Number of lines multilooked LLO Number of pixels multilooked 100 cox eoe eoe Fe Hee ee RA EERE ERER ie te o Fee AR AECA ce eR eoe eoe Ie AAA AA ie OK NC te AE Te oe e ANA te ede x End_interfero _NORMAL KEKKKKKKKKKKKKKKKKKKKK KKK KKK KKK KKK KK KKK KKK KKK KKK KKK KKKKKKKKKKKKKKKK The complex output data file can be viewed with e g in Matlab with a following script fid fopen Output cint raw r eime Erersacd Exe 100 22201 ilio S927 57 m fcloseifrd s zeeeuci Cilimes le2esilwe cum 1 15 Golsqoaicce emi 22822812 elm il 17 cint realpart ixcplxpart phase angle cint imagesc phase colombia Output may include the matrix with the reference phase For flatearth correction this normally resembles a plane and is not very useful output 75 pixel range 4 Buffer I Buffer II Rest Figure 23 2 Use of buffers in implementation of computation interferogram 23 3 Implementation The following is computed in buffers 1 Read in buffer of complex master and
138. k ck ck ok ck ck ok ok ck kk ck Sk Sk kk ko ko ko End comprefpha NORMAL Ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck Ck Ck ck Ck Ck Ck Ck Ck Ck CC C CC C CC CC CK CIC CK CCCII Sk Sk Sk S Sk S Sk Sk Sk KA kc kc kc kckok f 83 3 dE Up c th IRS o3 a A S3 99 CO S d amp 3 dc e Sy dis w RAS O IS XU RS SS SI quy jt OG do e e In the log file some statistics are given for the errors observation minus estimated value These errors should have a maximum of 0 1 phase cycle The can be plotted with GMT or someother package to evaluate the difference between the polynomial and the observations Here also the standard deviation per estimated coefficients is given This std seems to be too large but an error in the computations could not be found The polynomial can not easily be visualized at the moment It is normalized not evaluated in this step In the old method of the interfero step see Chapter 23 there was a card to output the reference phase polynomial because it was evaluated there anyway It seems logical to add a card for outputting the wrapped or unwrapped reference phase in the step subtrrefpha Chapter 25 which likely will be added in one of the cumming releases 79 Chapter 25 SUBTRREFPHA In this chapter the processing of step SUBTRREFPHA is described This step requires the steps INTERFERO and COMPREFPHA for obvious reasons In this step the reference phase of a mathematical body el
139. lipsoid is subtracted from the complex interferogram This is done by complex multiplication with the conjugated written symbolically as follows I I cos Rg isin Rg 25 1 Where is the complex interferogram denotes pointwize multiplication and Ro is the reference phase for a certain point 25 1 Input Cards SRP_METHOD polynomial exact method selector for subtraction of reference phase polynomial evaluates the polyno mial computed in step COMP REFPHA exact computes the reference phase explicitly for each pixel and subtracts it Computations are done by evaluation system of 3 equa tions foreach pixel SRP_OUT_CINT cint minrefpha raw filename of output datafile for complex interferogram of step subtrrefpha SRP_MULTILOOK 11 multilook factors in line azimuth and pixel range direction SRP DUMPREFPHA OFF ON This card specifies to dump the reference phase as a complex real4 file containing the evaluated reference phase polynomial just as it would have been subtracted from the complex interferogram multilooked The amplitude should be equal to one by definition WARNING the reference phase is not subtracted only dumped if this card is specified If you want to study the reference phase for different ellipsoids compile different versions of Doris changing the parameters in the file refsystem and use these executables to generate the reference phase SRP OUT REFPHA refphase raw name of output file r
140. ltering The area was relatively flat and the perpendicular baseline was approximately 175 meters 68 22 3 Implementation 22 3 1 Method porbits The frequency shift A f between the master and slave range data spectra equals c B _ tan A c Ad run UC NET cM Mcr MN 22 1 r tan a A tan 0 a A tan 0 a Where A0 0 05 a is the local terrain slope w r t the ellipsoid c is the speed of light 0 is the local incidence angle A is the radar wavelength r is the slant range ground to master The approximation is used in Doris Of course the sign of B1 or A0 is important to filter the correct side of the spectra Note that a 23 Af gt oo 22 2 The local incidence angle is computed with the dot product of vectors P and P M See also Gatelli et al 1994 The algorithm in Doris works as e While there is a line in the overlap get next line for master and slave e Get block of FFT LENGTH pixels e Compute viewing angle perpendicular baseline delta theta for middle pixel of block e Compute frequency shift by equation 22 1 and compose filter of rect and hamming e Filter master and slave e Write block back for last block only partially 22 3 2 Method adaptive After the resampling of the slave on the master grid is performed this algorithm can be used The local fringe frequency is estimated using peak analysis of the power of the spectrum of the complex interferogram The
141. ly a lot of detail is lost Figure 29 5 Magnitude of filtered complex inter ferogram Method spectral A pointwise mul tiplication in the spectral domain by a 32 point hamming filter was used a blocksize of 32 and an overlap of 4 96 Figure 29 4 Phase of filtered complex interfero gram Method spatialconv A spatial convolu tion with a kernel 1 4 9 4 1 was used Clearly a lot of detail is lost Figure 29 6 Phase of filtered complex interfero gram Method spectral A pointwise multiplica tion in the spectral domain by a 32 point hamming filter was used a blocksize of 32 and an overlap of 4 Figure 29 7 Magnitude of filtered complex inter Figure 29 8 Phase of filtered complex interfer ferogram Method goldstein parameters used ogram Method goldstein parameters used are alpha 0 5 smooth 3 overlap 4 This filter are alpha 0 5 smooth 3 overlap 4 This filter seems to preserve most detail seems to preserve most detail kernel 1 z H 29 1 This becomes 1 1 1 1 T 111 29 2 1 1 1 The blocksize for the convolution is chosen as high as possible A 2D kernel can be specified in an input file Only odd sized kernels can be used but simply add a zero to an odd kernel If a real4 matrix containing phase should be convoluted by a certain kernel first convert this real4 to a complex real4 matrix Do this either by computing the phase for complex umbers with amplitude 1 or by setting th
142. master and slave FC WINSIZE 32 32 44 Offset vectors Farga mg mn wp ze a um Figure 17 1 Plot produced by the command plotoffsets interferogram out 11 6000 21 1000 0 6 Out data 1393 raw keycard FC PLOT 0 6 BG The magnitude is plotted in the background Correlation is in dicated by the size of the circles estimates with a correlation below 0 6 are filtered out The size of the correlation window Recommended is 64 64 FC ACC 44 The search accuracy for the maximum correlation Adviced is 8 8 total search area is from Acc to Acc for FFT methods this must be a power of 2 In the logfile after step COARSECORR the variation of the initial offsets w r t the estimated values can be seen If this variation is larger than 1 1 is normal for ERS1 2 SLC images then one should select a bigger window and a larger search accuracy FC INITOFF 00 COARSECORR The initial offset between master and slave COARSECORR indicates that the results of the step COARSECORR are used FC OSFACTOR 16 The oversampling factor for the harmonic interpolation of the correlation Recommended is 32 to co register the images within a tenth of a pixel FC PLOT threshold 0 4 NOBG BG 45 Call gmt script plotoffset to plot results and to view with gv An example of a plotis given above This script gets the section with estimated fine offsets from the interferogram result file The argument threshold filters out estimates with a correlati
143. me The filename of the SLC data file This is the only file required for method ASAR EN VISAT S_IN_LEA filename The filename of the SLC leader file Not used for method ASAR ENVISAT S_IN_VOL filename The filename of the SLC volume file Not used for method ASAR ENVISAT and TSX TERRASAR X S IN NULL filename The filename of the SLC null file This may be a dummy name since it is not used 30 Chapter 10 S_PORBITS In this chapter the processing of step S PORBITS is described It is actually the same as step M PORBITS but then for the slave image See chapter 4 M PORBITS for more detailed information on this step 10 1 Input Cards S_ORBDIR directory name the tar archive directory name for the Delft Orbital Data Records S_ORB_INTERVAL 1 Time interval between data points Card M ORB INTERVAL has the same effect It is not possible to have a different S ORB INTERVAL S ORB EXTRATIME 3 Time before first line and after last for extra datapoints data points Card M_ORB_EXTRATIME has the same effect It is not possible to have a different S_ORB_EXTRATIME S_ORB_DUMP dt Dump interpolated t x y z to ascii file slaveorbit dat 31 Chapter 11 S CROP In this chapter the processing of step S CROP is described It is the same as step M CROP but then for the slave image See chapter 5 M CROP for more information on this step 11 1 Input Cards S CROP IN filename Filename of the SLC data file S CRO
144. memory considerations First complex interferogram is computed as in INTERFERO and the norms of the master and slave images are computed Then a shifting window of size COH WINSIZE is used to estimate the complex coherence see Annex D The coherence is computed with a function of the matrix class This function returns only the lines of the input which can be computed due to the edge of the the estimator window Then this is multilooked requires number of lines to be a multiple of the multilook factor Therefore the buffer should contain an overlap with the previous one 91 Chapter 29 FILTPHASE This chapter describes the processing of step FILTPHASE This step can be optionally used to filter the latest complex interferogram in order to reduce noise e g for visualization or aiding the phase unwrapping It is probably best run after step SUBTRREFPHA A lot of warnings can be generated if an image containing a lot of zeros is processed These warnings can be ignored The method goldstein is described in Goldstein and Werner 1998 Basically the fringes become sharper after filtering because the peak in the spectrum caused by the fringes is given a higher relative weight Method spatialconv simply is a spatial convolution with a certain kernel function e g a 3 point moving av erage Method spectral is a multiplication of the spectrum with the kernel specified in an input file e g a spectral low pass filter LPF 29 1 I
145. metry Shewchuk 2002 Shewchuk J R 2002 Delaunay refinement algorithms for triangular mesh generation Computational Geometry Theory and Applications 22 1 3 21 74 Touzi et al 1996 Touzi R Lopes A and Vachon P W 1996 Estimation of the coherence function for interferometric SAR applications In European Conference on Synthetic Aperture Radar K nigswinter Germany 26 28 March 1996 pages 241 244 Wessel and Smith 1998 Wessel P and Smith W H F 1998 New improved version of generic mapping tools released EOS Transactions AGU 79 47 579 Zebker et al 1994 Zebker H A Rosen P A Goldstein R M Gabriel A and Werner C L 1994 On the derivation of coseismic displacement fields using differential radar interferometry The Landers earthquake Journal of Geophysical Research 99 B10 19617 19634 121 Annex A What s new A 1 Version 4 02 The script doris process reset sh is renamed to doris rmstep sh and is updated for Mac OS X plotcpm is updated to plot offsets rather than the offset residuals Another script is made available under the name plotcpm residues to plot the offset residuals Bug fixes Step COMPREFDEM when multilooking applied during computation of reference DEM phase the out put file was smaller than the expected size Step COMPREFPHA fixed index overflow in sinus look up table A 2 Version 4 01 Handling of file sizes larger than 4GB on both 3
146. ming error based on correlation between mas ET molta dim amp S_CROP 11 See M_CROP SOVS 12 See M_OVS COARSEORB 13 Compute the translation between master and slave with the orbits precision 30 pixels COARSECORR 14 Compute the translation between master and slave on pixel level by correlation technique Spectral filter for master image in azimuth line direction S FILTAZI 16 Spectral filter for slave image in azimuth line direction FINE 17 Compute translation vectors over the total image on sub pixel level RELTIMING 18 Estimation of relative timing error between master and slave based on fine coregistration DEMASSIST 19 DEM assisted coregistration COREGPM 20 Compute the actual transformation 2d polynomial model for the alignment of the slave on the master image RESAMPLE 21 Resample the slave image according to the transformation model of the coregistration FILTRANGE 22 Spectral filter for master and slave image in range pixel di rection INTERFERO 23 Compute the complex interferogram COMPREFPHA 24 Compute the reference phase of the ellipsoid to be subtracted from the interferogram polynomial SUBTRREFPHA 25 Subtract the reference phase of the ellipsoid from the interfer E game ee ate ieee COMPREFDEM 26 Compute the reference phase of a DEM to be subtracted from NS e gram GEOCODE 33 Geocode the pixels convert pixels from the radar coordinate system to a ear
147. n We compute the synthetic amplitude for a given master acquisition using orbital information and topographic data SRTM can be used to obtain the topography in 3 arcseconds Near global 90 m at the equator or 1 arcsecond USA only resolution The input DEM file must have the byte order of your platform in order to extract correct elevation value see SAM IN FORMAT option for details The DEM should be in the WGS84 system same as the orbit ephemerides The Doris distribution contains the utility construct dem sh to download and prepare SRTM data see Section C 2 12 6 1 Input Cards SAM IN DEM filename filename of input DEM gtopo30 File is assumed to be stored in a raster Major row order from North to South line by line See also internet links at Doris home page for available DEMs SAM IN FORMAT 12 12 BIGENDIAN R4 R8 format of input DEM on file signed short for gtopo30 or real4 or real8 input ma trix is raw binary data w o header endianness of host platform is assumed except for I2 BIGENDIAN SAM IN SIZE 6000 4800 Number of rows and columns of input DEM file Default is set to tile w020n90 DEM SAM IN DELTA 0 00833333333333333333 deltalon Grid spacing of input DEM in decimal degrees latitude longitude Default is equal gridspacing default set to tile w020n90 DEM SAM IN UL 89 995833333333333 19 995833333333333333333 Coordinates of UL upperleft corner in decimal degrees latitude 90 90 longitude
148. n D 8 gives some information on the formats of the images D 1 Constants Constants used in the processing can be found in the source files constants h and refsystems h The main parameters in constants h are const real8 SOL 299792458 speed of light in m s const real8 ERS 1e 13 small number const int32 NaN 999 Not a Number const real8 BI 4 xatan 1 The main parameters in refsystems h are actually only WGS is used for now const real8 WGS84_A 6378137 0 semimajor axis wgs84 const real8 WGS84_B 6356752 3 semiminor axis wgs84 const real8 GRS80_A 6378137 0 semimajor axis grs80 const real8 GRS80_B 6356752 3 semiminor axis grs80 const real8 BESSEL_A 6377397 155 semimajor axis bessel const real8 BESSEL_B 6356078 963 semiminor axis bessel const real8 RADIUS 2 5 WGS84_A WGS84 B for sphere pol2car D 2 Baseline The basic configuration of InSAR is shown in figure D 1 There are different representations for the baseline see Figure D 2 Conversions between baseline representations The baseline parameters can be computed when the statevectors of the points M S and P master slave and point on surface are known The distance between the points x and y is denoted by d x y and the sharp 144 Figure D 1 Geometric configuration for InSAR R are the range vectors to the corresponding resolution element The statevector of the reference satellite is denoted by p h denotes
149. n square brackets The input file consists of a header and a tail In the header the general cards are placed see Chapter 2 and in the tail the cards specific to a certain step are placed described in the other chapters The order of the cards is not restricted except the STOP card though we advice to group them by processing step Blank lines are allowed in the input file but the line counter will not function properly in that case which does not affect the processing in any way We advice to place a comment on otherwise empty lines If accidently a certain card is used more than ones then a WARNING is generated and the latter one is ignored this behavior is not guaranteed not true for PROCESS cards The case of cards and parameters is not restricted We advice to use UPPER case except for comment or c for cards and lower case for parameters Text after the last expected parameter is ignored Be careful with putting comments in like this if the number of parameters may be varied for a certain card We normally do not make one big input file but we use several small ones for a group of processing steps See also the run file in Annex B The examples in the next chapters will make this more clear All keywords are described in this manual and also in the interactive helpdoris script and in the run script It can occur that some keywords are not mentioned in the latter two To obtain all possible cards give the command gre
150. ncy shifts larger than half the bandwidth A factor of 4 for example might give a better estimate since the interval between shifts that can be estimated is in that case halfed fixed FFTLENGTH RF_WEIGHTCORR ON OFF For method adaptive In peak estimation weight values to bias higher frequencies The reason for this card is that the low frequencies are for small OVERSAMPLE factors aliased after interferogram generation The deweighting is done by a dividing by a trian gle function convolution of 2 rect functions the shape of the range spectrum Effect of this card may be neglectable RF_OUT_MASTER master rfilter Output data file name of master RF_OUT_SLAVE slave rfilter Output data file name of slave RF OUT FORMAT cr4 ci2 Output data format for master and slave file Example input e e comment ADAPTIVE RANGE FILTERING Ci RF METHOD adaptive RF FFTLENGTH 128 f 2500 m RF_NLMEAN 15 odd RF_THRESHOLD 5 SNR RF_HAMMING 0 755 alpha RF OVERSAMPLE 4 RE WEIGHTCORR OFF RF OUT MASTER Outdata 3397 rfilter RF OUT SLAVE Outdata 23070 rfilter RE OUT FORMAT eq 67 20 Lig Zm 0 T T T C T T 0 0 01 02 03 04 05 06 07 08 09 1 0 E MT Nov 28 15 40 57 2000 Figure 22 1 Frequency histograms of the FINE coregistration correlation values At 301 locations of a image the fine coregistration
151. nd after last line to output ephemerides Since interpola tion is done with natural cubic splines it is advised to have at least 3 extra data points before the first and after the last line To use a single polynomial of degree 3 for in terpolation of the orbit for the full scene select a time interval of 20 seconds and for example extra time of 200 seconds M ORB DUMP delta t Write interpolated orbit to ascii output file masterorbit dat With delta t seconds interval between ephemerides Time interval between t 0 and t N of the precise ephemerides output is tx yz xdot ydot zdot xddot yddot zddot f compiled with DEBUG defined then also the matrices for spline interpolation are dumped Example of the cards for this step Ci e comment ERRORES e M ORBDIR data delftorbits ERS1 M ORB INTERVAL il M_ORB_EXTRATIME 6 c M_ORB_DUMP opa 4 2 Output Description If a normal termination of this step then the process flag at the start of the result file is switched to 1 precise_orbits ip The output of this step is written in the section precise datapoints This section looks like the following It is important that all lines are present following NUMBER OF DATAPOINTS 23 okbckck ckckockc kc ko RAE REX AA EA RE KK kok RARE kokok kok Kok kck kokck kc k kok kok AER ERA RAE RA k Kok x Start precise orbits kk ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck Ck Ck ck Ck Ck Ck Ck Ck Ck C C C C CC CC CC CK CIC RARA
152. neradlnputCards 2 000 eee eee eee 2 2 Example General Input Cards M_READFILES Sal JdInp tGardS s zou ae wee dS eS yin e Rm m Red Ud US s 3 2 Output Description uuo REPRE EE wi el S 3 3 Implementation koe aa RET RT 3 3 1 Changes for X86 platforms lens M PORBITS 41 ImputCards 0000000000000 08 4 2 OutputDescription sn i aa a a o 4 3 Implementation iuret encre oe a ee eee a ay Say RT TRE Sar Die 1s Fas M CROP 5 1 Inp t Cards z s s Ete HE RI RIDERE ASIE 5 2 Output Description a ecer e ae eb ee ere we ew we ee M_SIMAMP 6 4 INPULGATAS sree ceo ale Sa mm NA ee ee eee ee 6 2 Output Description ee 6 3 Implementation 00020000045 M_TIMING 1 Inp tGards sx ghee db bere Hh a EE RO Eos T2 OUIPUUDESCHIPUOM isai en rao SA ee se niti Boden coded eee dC Hus 7 3 Implementation llle M OVS 8 1 Inp t Qards opu 34 Oe P O IE ER eg 8 2 Output Description o ado x a a ae Gee os ms 8 3 Algorithm 2 6 9 0 65 e vod dr eI IP UYUE Y odes moz xod S READFILES 9 Inp t Cards s s sutor ulum s RERO SESS ada a 10 S PORBITS 10 1 Input Cards sass 9e sock ex R9 RH 93 3t we cm ow eee E 11 S CROP 1tj InputOards 2226999239 0b bd B ARA oh om m x xS 12 S OVS 12 1 InputGards 2 424640 eco os fm ee ee AAA n m 13 COARSEORB 19 1 Input Cards eave ea eho hbase eed aa
153. nj a2 ib2 EXAMPLE cpxdiv cint raw cint2 raw division raw 141 C 2 9 cpxconj Program to take conjugate of a complex file Program cpxconj take conjugate of a given complex float file USAGE cpxconj infile1 outfile infile 12 contain complex values a ib outfile contains conj a ib a ib EXAMPLE cpxconj cint raw cint raw conj C 2 10 floatmult This utility multiplies a complex float file by a scalar PROGRAM floatmult pixelwise float multiplication of a complex float complex file To be used to scale float files or magnitude of complex files SEE ALSO cpxmult flapjack cpxfiddle USAGE floatmult infile1 factor 2 EXAMPLE floatmult cint raw 1 47 C 2 11 wrap With this utility you can wrap your interferogram to arbitrary interval instead of pi pi Can be used for example to make fringes correspond to e g 1 cm displacement PROGRAM wrap wraps float binary file to interval a b USAGE wrap infile a b ofile EXAMPLE wrap interferogram raw 4pi 4pi interf4 raw default outputfile infile wrap default interval a b pi pi C 2 12 construct dem sh This utility downloads merges and fills voids of SRTM data based on coordinates of your area of interest Only basic Linux Unix commands and GMT are used so make sure you have GMT installed PROGRAM construct dem sh ownloads merges and fills voids of SRTM data based on coordinates of your area of interest USAGE c
154. not have a constant size with this method but varies between winsizeL P and 5winsizeL P 38 Chapter 15 M_FILTAZI In this chapter the processing of step M FILTAZI is described This optional step filters the spectrum of the master in azimuth direction The part of the spectrum that does not overlap with the spectrum of the slave is filtered out This non overlap is due to the selection of a Doppler centroid frequency in the SAR processing which normally is not equal for master and slave image This step can in general best be performed after the COARSE coregistration and before the FINE The coarse offset in pixel direction is used to evaluate the polynomial for the Doppler Centroid frequency The FINE steps can benefit a lot from this filtering TODO add plots By processing the RAW data to SLC at the mean Doppler centroid frequency this step can be avoided in the InSAR processing chain For ESA SLC images this cannot be done obviously Normally the step S FILTAZI is performed at the same time requires a PROCESS S FILTAZI card in the input file see chapter 16 However we kept this two seperate steps to be able to only filter the slave images in a large stack all slaves coregistered on the same master image This has the advantage that for each interferogram of the stack not a large file is created for the master The disadvantage of not filtering the master of course is that a small part of the spectrum of the master is not share
155. nput Cards PF METHOD goldstein Select goldstein method goldstein or spatial convolution spatialconv with cards PF KERNEL and PF IN KERNEL2D or spectral filter spectral with cards PF IN KERNEL2D PF BLOCKSIZE and PF OVERLAP For more info see implemen tation section PF_OUT_FILE cint alpha filtered Output filename for complex real4 file with filtered phase goldstein filter where alpha is substituted For method spatialconv default is cint filtered PF IN FILE filename numlines Optional filename of complex real4 inteferogram mph file to be filtered instead of the default which is obtained by reading the products result files Also specify the num ber of lines in this file as second argument This card is included to be able to filter files without having to create dummy result files to trick Doris For now the interferogram to be filtered must be complex real4 PF ALPHA 0 2 This card is for method goldstein only Alpha parameter for filtering This parameter must be between 0 no filtering and 1 most filtering The card PF_KERNEL influences this value since a higher smoothing relative decreases the peak and thus the effect of alpha PF OVERLAP 3 92 This card is for method goldstein and spectral only Half of the size of the overlap between consecutive blocks and buffers so that partially the same data is used for filtering The total overlap should be smaller than PF BLOCKSIZE If this parameter is s
156. nsequently the oversampling is performed line by line too In particular given an input signal Input signal xxxx and an oversampling ratio of let s say 3 I first generate a zero filled copy of the input signal zero filled Input signal x00x00x00x then I convolve the zero filled Input signal with the interpolation kernel obtaining the output signal CUASI GIAO RAPHP where s are the new samples Bert Kampes implemented the azimuth oversampling In azimuth a 6 point raised cosine kernel is used The kernel is normalized also normalized the range kernel typically 16 point sinc 29 Chapter 9 S_READFILES In this chapter the processing of step S READFILES is described It is the same as step M READFILES but then for the slave image See chapter 3 M READFILES for more information on this step 9 1 Input Cards S IN METHOD ERS ASAR ENVISAT RSAT RADARSAT ATLANTIS JERS ALOS TSX TERRASAR X Method selector to read ERS ENVISAT RADARSAT JERS ALOS or TERRASAR X header Note that both master and slave need to be acquired by the same sensor in principle JERS simply uses ERS programs ATLANTIS sar processor uses the ceos reader for RSAT and will write this in the Product Type Specifier field RSAT must be tested problems may be orbit data In later steps the Product field is read and the CROP step uses the appropriate function automatically Envisat ERS JERS RSAT ATLANTIS S_IN_DAT filena
157. nstallation of the TERRASAR X reader For the use of Terrasar X data the following additional packages are required on your system note the minimum version numbers e gdal version gt 1 44 See http gdal org for more information e python version gt 2 2 e libxml2 version gt 2 7 2 e python Ixml version gt 2 0 e libxslt version gt 1 1 15 In case your system does not meet these requirements and you cannot update you can try an alternative script In that case the requirements are e gdal version gt 1 44 See http gdal org for more information e python version gt 2 2 e libxml2 version gt 2 6 30 e python Ixml version gt 1 3 3 1 e libxslt not required To use the alternative script go to your Doris bin directory and do cp tsx dump header2doris no xpath py tsx dump header2doris py to overwrite the original script B 1 5 Starting Doris Now after a rehash we can run the Doris software and start InSAR processing The run script in the bin directory can be customized by setting the environment variables EDITOR and PAGER Note that it is highly adviced to install the utilities see annex C that can be found in the download area and GMT visualization To use the Delft precise orbits it is convenient to have getorb on your system but one could also use the online version from http www deos tudelft nl ers precorbs B 1 6 Installation of utility scripts In the Doris arc
158. ocation correct make a symbolic link as root 129 Table B 2 Troubleshooting 2 Problem Solution Run script doesn t work properly Run script editing option doesn t work properly not to my liking Run script viewing option doesn t work properly not to my liking Doris crashes at getting the Delft orbits Doris gets zero correlation in FINE coregis tration Doris cannot generate postscript while coregistration Script plotoffsets or plotcpm does not work Doris leaves temporary files in the working directory What is this file CPM data after FINE coreg istration Doris crashes on reading input Path not set to bindir try if command doris v returns ver sion number or GMT is not installed Or directory not writable chmod Use the environment variable EDITOR to specify your preferred editor for csh include in your startup file cshrc for example a line setenv EDITOR vi Use the environment variable PAGER to specify your preferred viewer for csh include in your startup file cshrc for example a line setenv EDITOR vi Is getorb installed Check location of orbdir Doris input card Obtain the system command Doris parses from the stdout INFO repeat it from the command line Is path to files correct in master res slave res Compile a debug version of doris and run that to find out what goes wrong Is GMT installed Check script from prompt by using command that is written as INFO
159. oducts result file is switched to 1 at successful exit fine_coreg i Example of output of this step products result file KEK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KK KKK ko ko ko KK _Start_fine_coreg KEKKKKKKKKKKKKKKKKKKKKKKKKKKK KKK KKK KKK KKK KKK KKK KKK KKKKKKKKKKKKKKEKKK Oversampling factor 16 Number_of_correlation_windows 100 Number posL posP offsetL offsetP correlation 0 2 2 241 06 3512 ORG il 29 2 3 J 241 12 306 Opis 2 20 447 241 12 Sl 0 44 SKIP SKIP 99 4733 761 241 06 3 44 OFS 100 4733 Qt DUANE 3 44 0 14 EXAXEXEXEXAXAENELKERRERE kckck ck ockokok kokcokok kokckok ERA RAE RRA ERA ARA AAA ERA kokck ck ko End fine coreg NORMAL KEKKKKKKKKKKKKKKKKKKKKKKKKKKKK ck ckckck ck ckckckckckckckckckckckckckckckckck ck ck ck ck ck ck ck ckckckck k ck kk In the logfile addition information is given 17 3 Implementation The current names for the master and slave image are read from the result files crop section Here also the dimensions of the files are read This can be checked with the debug version of Doris Doris can be tricked to coregister other complex files e g complex interferograms for 4 pass differential interferometry by substitution the right parameters in that section 17 3 1 magspace The computations are similar to the COARSECORR magspace method 17 3 2 oversample See source code 17 3 3 magfft The correlation is computed at pixel level similar to
160. of input DEM per buffer cut to the interfero gram window Info on these files is written as DEBUG CRD OUT DEMI filename Request optional debug output to float file of interpolated input DEM cut to the interfer ogram window Info on these files is written as DEBUG CRD OUT FILE refdem raw Filename of output radarcoded DEM CRD OUT DEM LP filename Request optional output of DEM in radar coordinates CRD OUT H2PH filename Request optional output of height to phase factors Most of these parameters can be found in the HDR file of gtopo30 DEM s Example input section REFERENCE DEM RD IN DE E final wanaF2835 dem RD IN FO AT r4 default is short integer RD_IN_SIZE 3601 3601 RODEN DEA OR OW CS55 5 Sm 01 0 0 008 9 3 9 3 RD_IN_UL AQ 27 RD IN NODATA 32768 RD OUT DEM LP 42408 22735 demlp BDOOUTCRILE 42408 22735 demphase RD OUT H2PH 42408 22735 h2ph CPG nme mene me 26 2 Output Description At successful exit the process control flag is switched on comp_refdem il 84 Figure 26 1 Interferogram of radarcoded DEM for area described in section 23 0 4 The output in the products result file looks like KKK ck ck ck Ck ck KK KKK KKK KK KKK KKK KKK KKK KKK KKK KKK AAA AAA RARA ARA KK KKK KK Start comp refdem kc ck ck ck Ck ck kk ck Ck ck ck Ck ck kc ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck koc ck ck Ck ck ck AAA RR
161. oks Snaphu option c Refer to snaphu manual for more information UW_SNAPHU_LOG filename Output log file name for snaphu option l Refer to snaphu manual for more information UW SNAPHU INIT MST MCF Output log file name for snaphu option l Refer to snaphu manual for more information UW SNAPHU VERBOSE ON OFF snaphu option v Refer to snaphu manual for more information 30 2 Output Description At successful exit the process control flag is switched unwrap il The section for the unwrapping in the result file for the interferogram looks like the following for method TREEF file name and format are used later kk Ck ck ck ck ck ck ck ck ck ck ck ck ck KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK KKK Sk KKK Sk S ke Sk Sk kc ko KKK AAA x Start unwrap Ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck Ck ck ck Ck Ck Ck Ck Ck C C C C CC C C CC CC CIC AAA RRA RARA Sk Sk S ke S ke ke KK kc k kc kc KKK Data output file Data output format Data output file regions Data output format Last line w r t original master Multilookfactor azimuth direction Multilookfactor range direction Program for unwrapping Output program for unwrapping Delta lines for seed Delta pixels for seed Number of patches used First line w r t original master Sci o ali clas eta Last pixel w r t original master Outdata uint raw real4 Outdata regions raw short int 2B AO ORE ZOS SOn 700 1L 0 2
162. olution cells As a result a non square correlation window may be more suitable For example for ERS1 2 we apply a 256 lines by 128 pixels correlation window This step updates the master acquisition azimuth and range times respectively Therefore it affects steps whenever the master timing is used such as the coregistration using a DEM DEMASSIST step and the computation of reference phases COMPREFPHA and COMPREFDEM steps 7 1 Input Cards MTE_METHOD magfft magspace Method selector for this step Either perform the correlation computation on the magni tude images in the space or in the spectral domain MTE IN POS filename Input filename for ASCII file with positions in original master system to place windows for correlation computations MTE NWIN 16 Number of windows to be distributed over the total image to estimate the offset Should be at least 5 or so because the most consistent estimate is selected This card is ignored if MTE IN POS is set Only 1 large window could be used e g of size 1024x1024 MTE_WINSIZE 256 128 Size of the window in lines pixels For method in space domain it defaults to 256 128 For method in space domain it is converted to odd numbers if necessary MTE_ACC 32 32 25 ONLY for method in space domain Accuracy to search within for maximum correlation For fft method it automatically equals half of the MTE_WINSIZE In case of magspace the DEM window size is extended by 2xMTE_ACC MTE_INITOFF
163. omial in the header file 63 For a correct interpolation either the spectrum of the data has to be shifted to zero or the interpolation kernel We shift the kernel in azimuth In range the spectrum is centered for satellite data A simple derivation shows how the kernel should be shifted Suppose we have a signal s x exp ix2pixx xfdc prf and we want to interpolate this signal at xi 5 1 with a triangular kernel k x 0 such that the interpolated signal S xb dE al This means with x_0 0 1 0 9 and x_1 5 6 as implemented in Doris the kernel is formed as k x_0 triangle x_0 0 9 0 1 and we shift this with the MINUS sign by multiplication with t x_0 exp i 2pxx_Oxfdc prf Then the interpolated value at 5 1 equals sS r 5 1 le 00 53D see 50 1 m 5 a 0 9 eic 001 9 s 6 0 9xexp ix2pix 0 1xfdc prf xexp i 2pi 5 xfdc prf 0 1xexp ix2pix0 9 fdc prf xexp i x2pix6 xfdc prf exp i 2pix5 1 fdc prf 5 1 perfect interpolation On the contrary when we use the PLUS it follows that S_al 30 0 9x exp 1 2pix O0 1 fdc prf exp ix 2pix5 fdc prf 0 1xexp ix2pix0 9 xfdc prf xexp ix2pi 6 fdc prf 0 9xexp i 2pi x4 9xfdc prf O l exp i 2pix6 9 fdc prf NE s 5 1 wrong sign used See also Hanssen and Bamler 1999 21 3 1 Output formats Computations are done in complex float Casting these values to complex short format introduces an error O
164. on below this threshold The second argument BG or NOBG selects a call to cpxfiddle to generate a magnitude background while BG does call cpxfiddle See the script plotoffsets and the c program cpxfiddle for more information cpxfiddle can be downloaded from Doris internet pages The command is echoed to stdout as INFO which can be repeated outside Doris Before running the step COREGPM to estimate a transformation model it is very convenient to view a number of offset vectors above a correlation threshold to select the appropriate value for the card CPM_THRESHOLD Actually a background call is made to the script plotoffsets something like plotoffsets interferogram out 11 6000 21 1000 0 6 This command can be given from the prompt as well for different values of the threshold With a command like awk BEGIN for 1 100 i lt 25200 1 1 500 for 3 750 3 lt 5400 3 3 200 Pr Fa Vi a ah eeu the file for FC IN POS can be easily generated for a grid of locations Example input cards for this step e e comment _ FINE COREGISTRATION____ s FC METHOD oversample c FC METHOD magfft ul c FC METHOD magspace FC NWIN 0 1 number of windows FC_WINSIZE 64 64 size of windows HERACE 8 8 search window 2 n FC INITOFFE coarsecorr use result of coarse to compute first FC_OSFACTOR 32 oversampling factor 17 2 Output Description The process control flag at the start of the pr
165. on_pixels 3 KARA ok oe ode KKK eoe e ee oe oe ee ee KR he te ae lee ee ae eke e ee He ee Te KK He ER o He teak e te ee xx End_coarse_orbits _NORMAL DEI Te II Hee IEE ee I ee eo ooo Te ie RAS AR See eoe eoe eae RR RRA Ie ie ies He Nee Ue he x aie he x In the logfile some extra information is given such as the number of iterations The baseline parameters are not used but given here to make it possible to write scripts that grep these values if a lot of interferograms are processed of the same scene 13 3 Implementation The algorithm described in Annex D is used for the conversion between line pixel coordinates to the cor responding point P on an ellipsoid The Doppler range and ellipsoid equation This step consists of three steps basically 1 For the center line pixel of the master image compute the position x y z in system of the orbits of the point P on an ellipsoid 2 Based on the Doppler equation compute the position of the slave satellite corresponding to the point P on an ellipsoid and compute the line pixel coordinates in the slave system 3 The difference slave master between the line pixel coordinates is defined as the offset 35 Chapter 14 COARSECORR In this chapter the processing of step COARSECORR is described The offset in line azimuth and pixel range direction between master and slave is computed with an accuracy of about 1 pixel 1 offset for whole image The magnitude images a
166. onstruct dem sh project WE S N SRIM 1 3 lt ftp1 gt lt ftp2 gt lt ftp_user gt lt ftp_pass gt EXAMPLE construct dem sh netherlands 3 3 7 3 50 7 53 7 SRTM3 OUTPUT DEM final PROJECT dem preview srtm PROJECT ps Doris input input doris PROJECT 142 C 2 13 doris process reset sh This utility can be used to reset and clean up the processing entries in Doris result files such as master res slave res or interferogram res in order to allow re processing of a step or multiple steps Basicly it deletes any entry starting from the indicated step till to the very end of the result file and updates process control switches USAGE doris process reset sh doris process name lt res gt EX doris process reset sh coarse correl master slave res doris process reset sh resample slave res Doris process names by order 1 crop 13 resample 2 sim_amplitude 14 interfero 3 master_timing 15 comp refphase 4 filt_azi 16 subtr_refphase 5 filt range 17 comp refdem 6 oversample 18 subtr_refdem 7 coarse orbits 19 coherence 8 coarse correl 20 filtphase 9 fine coreg 21 unwrap 10 timing error 22 slant2h 11 dem assist 23 geocoding 12 comp coregpm 24 dinsar 25 extra C 3 Completes for tcsh users Complete commands are used in tcsh shell to complete commands by pressing the TAB key These completes can be added to the ones you already have Simply put them in your resource file likely cshrc or tcshrc
167. oris on a SGI platform Thanks to Kamini Kanta Mohanty From Kamini Kanta Mohanty lt mohantykk yahoo com gt Subject My Experience on DORIS To Bert Kampes lt kampes geo tudelft nl Cc doris usersOtudelft nl kkm_10 hotmail com Doris Listserver Delft Object oriented Radar Interferometric software From K K Mohanty Marine and Water Resources Division Space Applications Centre ISRO Ahmedabad 380 053 INDIA mohantykk yahoo com 7 July 2000 To Dear Doris Users At the outset would like to thank Mr Bert Kampes DEOS Delft University to make DORIS openly available 131 for download have downloaded DORIS 2 3 software sometime in the middle of May 2000 Subsequently have installed the same in SGI Silicon Graphics Octane w s with IRIX 6 4 o s have executed many steps in the s w not all would like to share my experience of installing s w in SGI machine Also got some queries May be a few suggestions which can be taken care of in future release am new to interferometry am based at Space Application Centre ISRO Ahmedabad INDIA am also an ITC Netherlands alumni Experience During Installation 1 complex type in SGI is a class not template supporting only double type had to separately bring in the complext h from GNU It works fine 2 Equivalents of all other includes are available but h extension has to be added for example
168. ormed to reduce noise Usually a ratio of line pixel 5 1 between the factors is chosen to obtain more or less square pixels 20x20m2 for factors 5 and 1 The resolution decreases of course if multilooking is applied Data buffering is applied for the memory considerations 28 1 Input Cards COH_METHOD include refdem refphase only 89 Figure 28 1 Coherence estimate using the REF Figure 28 2 Coherence estimate using the IN PHASE_ONLY old method The colorscale CLUDE_REFDEM new method The colorscale ranges from 0 to 1 ranges from 0 to 1 Method selector for coherence map generation INCLUDE_REFDEM computes the co herence map with subtraction of both reference phase and radarcoded DEM phase with the assumption that the last one has not been multilooked By selecting the REF PHASE ONLY method the coherence map is computed with subtraction of the only reference phase if available COH_OUT_CCOH filename Filename of output datafile for complex coherence image of step coherence one of COH_OUT_CCOH and _COH is mandatory COH_OUT_COH filename Filename of output datafile for real coherence image of step coherence One of COH_OUT_CCOH and _COH is mandatory COH MULTILOOK 102 Multilookfactor if no multilooking is required set this to 1 1 COH WINSIZE 102 Window size of shifting window for coherence estimation Example input section for this step E E comment product generation COH METHOD
169. orrelation_translation_pixels 3 KEKKKKKKKKKKKKKKKKKKKKKKKK KKK KKK KKK KKK KKK kckckckckckckckckckck ck ck ck ck ck ck ck AAA ck kk End coarse correlation NORMAL KKKKKKKKKKKKKKKKKKKKKKKK ckckckckckckockckckck ck ckckckckckckckckckckckckckckckckck ck ck ck ck k ck ck ckckckck k ck kk In the logfile the estimated offset is given for all windows 14 3 Implementation 14 3 1 Method magspace The implementation in the space domain requires an odd window size which is automatically forced not strictly necessary but this made the implementation a bit easier because the center of the shifting window is defined at a pixel For each location the zero meaned slave magnitude window is shifted over the zero mean master window and the correlation is computed see equation D 24 by computing all pointwize prod ucts and dividing by the norms of the particular windows 37 14 3 2 Method magfft The implementation in the frequency domain is more or less the same as in the space domain We only use FFT s to compute the products for the correlation see equation D 24 in an efficient way due to the fact that a convolution in the space domain corresponds to a multiplication in the frequency domain Input are the zero mean magnitude images The cross products are obtained by computing the pointwize product of the zeropadded master x conj slave A block function is used to compute the norms Note that the correlation window the overlap does
170. ove WGS84 This can be used in future to correct the unwrap ping for the integration constant Now it is only used if GEO card is used for cropping TIEPOINT lat lon hei Coordinates of a point in lat lon hei in WGS84 For now only informational The point is converted to pixel line coordinates and the interferometric phase is computed etc M RG T ERROR 0 0 Range timing error for master One way in seconds Use this card for example to cali brate the geo referencing using a corner reflector with known coordinates Can also be used to shift the DEM with respect to the interferogram in step COMPREFDEM A shift of one non mulitlooked or oversampled pixel corresponds to a one way timing error of 1 2 RSR For ers this is approximately M_RG_T_ERROR 0 00000002637 seconds By multiplication of the signal velocity speed of light 368 this amount in seconds can be converted to the slant range resolution i e pixel posting of 7 9 meter M AZ T ERROR 0 0 Azimuth timing error for master Use this card to account for timing errors in azimuth direction Card can be used to shift a DEM in azimuth direction Note that such a shift may indicate incorrectly estimated Doppler S RG T ERROR 0 0 Range timing error for slave See M RG T ERROR Since the geometry of the interfer ogram is related to the master this card has not a large effect S AZ T ERROR 0 0 Azimuth timing error for slave See M AZ T ERROR Since the geometry of the
171. p SUBTRREFPHA and sampled on the same grid see step RESAMPLE The files must have the same multilook factors and the same dimensions i e overlap exactly The perpendicular baseline of the topo pair should be larger than that of the defo pair to prevent that noise is blown up but this cannot always be controlled of course This step is performed in the defo pair processing tree First create a directory to run the topo pair processing until a unwrapped interferogram is obtained keep the master slave and products result files Then perform the defo pair processing After interferogram generation and flatearth subtraction start this step DINSAR specifying the location of the result files of the topo pair processing with input cards For 3 pass use a common master For 4 pass coregister the complex interferogram on the complex interferogram of the defo pair and then unwrap or first coregister master and slave on the master of the deformation pair To geocode the differential phase values geocode the topo interferogram and use the latitude longitude ma trices for the differential grid 31 1 Input Cards DIL OUT FILE differentialinterf raw Output filename for complex real4 file with differential phase in slant range system DI IN TOPOMASTER same as master result file card Specify this card if 4 pass differential interferometry is required Do not use this card for 3 pass or use the same name as the master result file for the de
172. p keyword readinput cc grep else 1 2 2 Outputfiles There are three ascii output files one for results of processing steps specific to the master one for the slave and one for the rest of the processing the products These files are referred to as master slave and product result file parameter files For example the wavelength of the sensor and the filename of the master image can be found in the master result file and for the slave parameters in the slave result file while coregistration parameters which aren t unique to a particular image can be found in the products result file These output files serve as input for Doris for running later steps Of course a step also can generate binary data output This is described in the following chapters The result files all consist of a header and a tail which grows with the processing In the header some general information and an overview of the processing is given with process control flags These flags do not imply a certain order By convention each processing step can be run only once 0 or 1 in the process control flag to avoid confusion on the correct latest results are This implies that a result section in the tail has to be deleted and the process control flag reset to 0 before running a step a second time In the growing tail the results of the processing is stored The result files are read again and the read parameters are used in the further processing In
173. pler The point P at the surface lies perpendicular to the orbit due to zero Doppler processing otherwise this equation has to be adapted with a slant angle 2 Range The geometrical distance to P on the surface is equal to the speed of light times the range time 3 Ellipsoid Force the point to lie on an ellipsoid The equations for the point P on the ellipsoid and the satellite S in its orbit are where x denotes x y z dr T Ts D 44 Ei s dv 0 D 45 Ez dz dx vjighttrange 0 D 46 2 2 2 Bia tuat eqep D 47 a a b To compute the coordinates of a point P on the ellipsoid corresponding with line and pixel p in the master image the following has to be done First the position of the satellite has to be computed assumed exact based on the line number and PRF I to azimuth time to interpolated position and the velocities for this time by interpolation Also the range time corresponding to the pixel number is computed based on RSR assumed exact Next the set of equations is used to solve for P x y z This is done iteratively by linearization which requires the derivative of the equations to x and approximate values for the unknowns the coordinates of the center A given in the SLC leader file converted to xyz on a sphere dz T D 48 151 ip E Ys d dz dE Ea _ i 5 A 2dx 2dy 2dz D 49 Ea Es 9E3 2x 2y 2z x dy z a a oF Solving this exactly determ
174. r azimuth direction 54 Range_direction Correlation 0 0 0 1 0 2 0 3 0 4 0 5 0 6 0 7 0 8 0 9 abs e abs e 0 0 0 1 0 2 0 3 0 4 0 5 0 6 0 7 0 8 0 9 Correlation Figure 20 3 Plot produced by plotcpm for the first run The absolute error estimated offsets minus observed offsets are plotted for the range direction Location windows wtests Range 399 798 1197 1596 1995 2394 2793 3192 3591 3990 14990 es Hae pb L i 1 14990 T7 ATi 480 481 lt 483 M67 Heo 472 4M54 M55 MSG H57 4 58 HMSO 13491 M42 448 A44 dS dG HAT M48 L13491 4434 436 20 H21 4422 423 ye 425 H26 407 4 409 40 ett m um 898 4401 03 11992 gg 4B85 886 4 B87 4 B88 en m 892 11992 B74 44876 4 877 4878 79 lt 880 Bee 6 B69 B49 B52 854 4857 10493 4839 pao V Seas 846 10493 x827 4 830 4 33 pu Are 18 RS 2805 B06 4807 8994 4 91 44903 494 4 897 00 H 8994 c 2281 1283 284 2286 p89 c 4MP69 270 deni 472 473 74 2256 0257 59 p60 p6i oper lt p65 5 E 74954 eit R49 252 4 253 L 7495 E 33 9235 f 41 PAZ N 4221 22 e24 25 26 31 N R10 4815 p18 xt 498 POO 401 202 203 205 og 5996 4 89 H90 lt H91 92 4493 4 97 5996 at 75 ATT a78 4664 4 466 44967 AO yn 474 957 4 462 4497 40 has 150 gt 4497 4428 4 3a 439 47 H18 c9 4423 ci 07 Sn cu2 efits iss 2998 94 295 4 B6 47 488 02 2998 82 83 B7 88 ea _ 1 3 eri 2 80 59
175. r module to the big selecting switch in main 3 Implement your module let in result file the output section end with same string as the other modules END NORMAL 4 Documentation author date description for users and code developers 5 Email the description and the total source to the owner of the mailinglist doris users amp tudelft nl And if approved we will include your functions in the next version of Doris F1 Formats The example source code explains which rules for commenting generally follow 0 DR RRR OK RRR GOR RK kok ek 1 ts16 2 3 truncated sinc 16 points 4 5 input 6 x X axis 7 Output 8 x y f x function evaluated at x 9 10 Bert Kampes 16 Mar 1999 11 SOOO ke ke oda ok e oe oe ke ad o ae le ale al GOI GG k k k I k k k k k kak ak kak o 12 matrix lt real4 gt ts16 13 const matrix lt real4 gt amp x 14 1 15 ifdef _DEBUG 16 DEBUG ts16 17 if x pixels 1 158 18 ERROR ts16 standing vectors only 19 endif 20 matrix lt real4 gt y x lines 1 21 for register int32 0 i lt y lines i 22 y i 0 sinc x i 0 rect x i 0 16 23 return y 24 END ts16 e start routine a block with date author description input output e end routine met END routinename e no block comments inside the function only things like Comment on block Mo Comment on something smaller ______ Indenting
176. re used correlation is computed in the space or spectral domain At a number of positions geometrically distributed or at positions read from an input file in the image the correlation between master and slave is computed for different offsets The offset with the highest correlation is the estimate for that position The approximate offset between the two images is set to the offset that most occurred over the positions so the one that is most likely Sometimes an estimated offset is totally unreliable for example for a position in a sea but the correlation is not very small The estimated correlation at a position is likely to be biased Therefor it would not be wise to use the offset between the two images based on the highest correlation values only but we use this consistency test instead 14 4 Input Cards CC_METHOD magfft magspace Method selector for this step Either perform the correlation computation on the magni tude images in he space or in the spectral domain CC IN POS filename Input filename for ascii file with positions in original master system to place windows for correlation computations CC NWIN 11 Number of windows to be distributed over the total image to estimate the offset Should be at least 5 or so because the most consistent estimate is selected This card is ignored if CC IN POS is set Only 1 large window could be used e g of size 1024x1024 CC WINSIZE 64 64 Size of the window in lines pixel
177. ri cos 32 19 Then compute cosine of angle mu across r1 in same triangle r pi p 2rhoip cos u 32 20 Unclear how to compute H exact for now use approximation set radius of earth at location of satellite equal to radius at location of P compute satellite height by Bowring s method xyz2ell then radius R of earth at phi lambda to satellite R pi Hest 32 21 Approximate H in this way H p Reosp 32 22 Compute error of this approximation preliminary study This will cause a bias and some trend in the height Because there likely already is a trend due to orbit errors this is not as bad as it might seem By using tie points a good height may be computed For now we did not implement a routine that uses tie points New way of computing H NOT implemented to difficult 1 compute in new system x y coordinates of 1 0 rho1 P 2 ellips equation in tje same system rotated over co latitude 3 snijpunt Bellips R 4 H rho1 Rq Problems with this method how do you know orientation of theta rotation of ellips to new system A few more notes A By 7 or 32 23 T So the reference phase has to be added again in order to compute Bpar otherwise Bpar is 0 001 m or so Processing 1 per line compute B alpha 112 2 per pixel a phi to Bpar b r known c compute theta exact d compute p e compute mu f compute H g compute h The idea is to comput
178. ript that generates the Makefile again Consider adding a library routine for FFT to doris since the internal one is not optimized speedwise Compilation fails due to int16 problem with Change code in ioroutines c routine checkrequest from var arg in ioroutines c int16 to int e g int16 N va arg arglist int16 gt int N va arg arglist int This happened with redhat 7 0 but not sure if Doris still functions completely correct It should have no influence on the computations only on user friend lyness i e warnings if Doris thinks certain steps should not be run Compilation fails because compiler cannot Change include files in source e g include lt cctype gt find include files to include lt ctype h gt Particularly older compilers do not know of the new standard e g cctype Compilation fails due to strptime Redhat Remove comment before DEF8 in the Makefile i e DEF8 7 0 D_NO_STRPTIME Remove object files with command make clean not required Compile again with command make Run script doesn t work at all Check first line of run script There the interpreter an ex ecutable program is specified for example bin ksh Does this program exist test with command which ksh For redhat 7 0 who doesn t include it in the distribution for some strange reason search the internet for a public domain ksh korn shell search at www google com on pdksh download Is the l
179. ry part Table D 2 The way complex files SLC data resampled slave complex interferogram are stored on disk Also refered to as mph format magnitude phase This is a major row order stored pixel interleaved file with 2 float 4B canals real imag 1st pixel 2nd pixel pixel P Tatline Real Imag Real Imag Real imag 2ndline RealImag Real imag __ Sine Real mag Realmag p TENER a Er pec inel Real imag Real mag Realimag After unwrapping of the phase the result can no longer be stored as a complex value because a complex number only can distinguish between phase values in the principal interval 7 Therefor a new format is used will be used Either the unwrapped phase is simply stored in a 4B float file similar to table D 2 without the imaginary part as are other files like the phi lambda and height matrices after geocoding or a hgt file is generated see table D 3 Particularly after unwrapping the conventions we use for this file are as follows The amplitude equal to that of complex intererogram and the unwrapped phase is stored for each pixel If there is no unwrapped phase the wrapped phase is stored and the amplitude is set to O for that pixel would prefer setting the phase to 0 and keeping the amplitude but we selected this format to keep in line with other software The amplitude is stored while it does not change after unwrapping to keep all information
180. s For method in space domain it defaults to 64 64 For method in space domain it is converted to odd numbers if necessary CC ACC 328 ONLY for method in space domain Accuracy to search within for maximum correlation For fft method it automatically equals half of the CC WINSIZE CC INITOFF 0 orbit 36 Initial offset for coarse co registration if the word orbit then the estimate of the step COARSEORB are read from the products result file and used If there are 2 numbers then these are used Example input cards for this step El s comment COARSE CORR COREGISTRATION s CC METHOD magfft default CMC CETT ET I9 magspace Hf ae Mec Ls c CC ACC 30 30 only for magspace CC_NWIN 21 number of windows CC_WINSIZE 1024 512 size of windows CC INITOFF orbit Vi tes result Oi exclus Ll for initial offset eI UNTER ERR use this if no precise orbits 14 2 Output Description At successful exit the process control flag is switched to 1 in the products result file If this file does not exist it is created ILRESFILE card coarse_correl al The output in the products result file resembles KEKKKKKKKKKKKKKKKKKKKKKK KK KK KK KKK KKK KK ARA KKK KKK KKK KKKKKKKKKKKKKKKKK _Start_coarse_correlation KEKKKKKKKKKKKKKKKKKKKKKKKKK KKK KKK KKK KKK KKK KKK KK KKK KKKKKKKKKKKKKKKKK Estimated translation slave w r t master Coarse_correlation_translation_lines eal Coarse_c
181. s seems to yield a more or less scaled version of the height of the ambiguity method The schwabisch being higher After rescaling with a factor 0 86 of the heights obtained by schwabisch method to the level of the ambiguity method the differences between both methods were only a few meters 116 Chapter 33 GEOCODE In this chapter the processing of step GEOCODE is described In this step the radar coded heights are converted to geocoded coordinates i e to a known reference system Input is the height file from the step SLANT2H Output are two files containing the latitude phi and longitude lambda corresponding to the height file These files can be further processed with programs cpxfiddle proj and GMT to create a DEM in a regular grid in any projection desired UTM for example See the bin directory and the shell scripts there for examples how to do this Likely for each specific application these scripts are best copied and adapted to your needs The interpolated matrices from GMT are in grd format that can be handled by Matlab etc 33 1 Input Cards GEO OUT PHI geo_phi raw Output file name for latitude GEO_OUT_LAM geo_lam raw Output file name for longitude Example input E comment GEOCODING e GEO_OUT_LAM Outdata lam raw GEO_OUT_PHI Outdata phi raw If you want to obtain the latitude longitude of the pixels in an interferogram that was created but you do not have a DEM in radarcoor
182. s and run script where a sh implementation is used instead of ksh linked csh to tcsh in the bin directory of my cygwin installation which also works fine and prevents changing the scripts i e n s bin tcsh bin csh Install GMT tools Install XFree86 and ghostview xv etc e i didn t install getorb and that may partly use a fortran compiler If someone has done this please let me know e For convenience make a symbolic link to your cdrom drive e g n s d cdrom You can then refer to SLC files on cdrom in the Doris input with cdrom SCENE1 B 6 List of files in archive The following directories are created after tar xvf Dorisv1 0 tar e Bin csh scripts helper programs e Source new Source code Doris For example the following files are in the archive Dorisv2 5 tar 11 July 2000 yee 5767 Jul 7 17 47 2000 Source_new makefile 413 22 8998 Jul 17 47 2000 Source_new constants h 413 22 19441 Jul 17 47 2000 Source_new conversion c 413 22 5641 Jul 17 47 2000 Source_new conversion h 413 22 127641 Jul 17 47 2000 Source new coregistration c 413 22 5020 Jul 17 47 2000 Source new coregistration h 413 22 20883 Jul 17 47 2000 Source_new filtering c NNNNNN 133 ANoOahWNMHO ak ks nk M Oo 111 413 22 1410 Jul 17 47 2000 Source_new filtering h r r 413 22 60217 Jul 17 47 2000 Source_new geocode c r r 413 22 3549 Jul 17 47 2000 Source new geocode h r r 413 22 1
183. s output is appended to the master result file KC ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck Ck Ck ck ck Ck Ck Ck ck Ck kk ck Ck Ck Ck C CK CIC CK CCCII Sk Sk e S ke Sk ke kk kc k kc kc kc kckok x Start readfiles KC Ck ck ck ck ck ck ck ck ck Ck ck kk ck ck Ck Ck Ck Ck C CC CC CC CC C CC CC CC CK CIC RARA RRA RARA Sk Sk S Sk RARA k ko ck ck kc kc kckok Volume file cdrom scenel vdf dat 001 Volume ID il Volume_identifier 0004093800014027 Volume_set_identifier 19950830 GAG Sei Check Number of records in ref file 26558 Product type specifier PIRODUIC W SIRS Il SYR SILC Location and date time of product creation IPAF 24 07 1998 Scene identification ORBITZ ISO DALE 30208795 Scene location FRAME 2781 LAT 40 94 LON 14 03 Leader file cdrom scenel lea 01 001 Scene centre latitude 40 9380000 Scene centre longitude 14 0270000 Radar wavelength m 0 0566660 First pixel azimuth time UTC SUSAUG TOPS OS 9I A59 Pulse_Repetition_Frequency actual Hz IETS 2020000 Total_azimuth_band_width Hz 1378 0000000 Weighting_azimuth HAMMING Xtrack_f_DC_constant Hz early edge 437 9780000 Xtrack_f_DC_linear Hz s early edge 7154 0000000 Xtrack f DC quadratic Hz s s early edge 380000000 00000 Range time to first pixel 2way ms 5 5458330 Range sampling rate leaderfile MHz 18 9624680 Total range band width MHz 15 5500000 Weighting range HAMMING Datafile cdrom scenel dat 01
184. sed filter for slave red is composed r r r r r 2E r r r r r 0 f 0 f 800 600 400 200 0 200 400 600 800 800 600 400 200 0 200 400 600 800 filtered spectrum for master filtered spectrum for slave 450007 T T z T T 7 15000 10000 109001 5000F 5000F l VARIA Emme J A ATO ne Bau AAA 800 600 400 200 0 200 400 600 800 800 600 400 200 0 200 400 600 800 Frequency Hz Frequency Hz Figure 15 1 Azimuth filtering for a master left and slave right SLC image frame 2781 orbit 1393 master ERS2 27 JUL 1995 and orbit 21066 slave ERS1 26 JUL 1995 The Doppler centroid frequency for the master is fp m 117 Hz constant for all columns for the slave f pc 425 Hz obtained from the result file read from SLC leader The mean Doppler centroid equals fpc 271 Hz Doppler centroid are indicated by dashed magenta lines x axis are frequencies from PRF 2 PRF 2 The azimuth spectrum was weighted with a Hamming window a 0 75 Pictures on first row original spectra for range column 101 red dashed line is a 51 point moving average The filtering middle row first de weights by inverse Hamming centered at the image Doppler centroid and bandlimited to the total azimuth bandwidth ABW 1378 Hz Next a new Hamming filter is applied centered at the mean Doppler centroid and bandlimited to ABW 2 focm foc 1070 Hz Obviously
185. sed by errors in the software or in the documentation 4 Users are very welcome to extend the capabilities of the Doris software by implementing new algorithms or improving the existing ones It is intended that if new software is developed based on Doris that this also is made available for free to the other users through us 5 We would appreciate if any addition or modification of the software would be announced first to us so that it can be included in the official next version of the software 6 Publications that contain results produced by the Doris software should contain an acknowledgment For example The interferometric processing was performed using the freely available Doris software package developed by the Delft Institute of Earth Observation and Space Systems DEOS Delft Uni versity of Technology or include a reference to Bert M Kampes Ramon F Hanssen and Zbigniew Perski Radar interferometry with public domain tools In Third International Workshop on ERS SAR Interferometry FRINGEOS Frascati Italy 1 5 Dec 2003 page 6 pp 2003 This document is typeset in ATEX2e Delft December 2008 Contents Preface 1 Introduction 1 1 Overview of the INSAR Processing n 1 18 Processingoorder 2 nm oL3 Aa 1 2 General considerations and conventions ee EAS MEI PTT CT 12 2 OutputflleS or aes ae beeeey dro Ry eee X EE ee Rs 1 3 Outlineofthisdocument llle General Cards 2 1 Ge
186. ses height ambiguity to compute the height schwabisch method uses polynomials to compare the phase with the reference phase rodriguez method uses geometry contains an approximation it is not clear how to compute a certain parameter S2H OUT HEI hei raw Output file name for computed height values S2H OUT PHI phi raw Only for ambiguity method Output file name for computed phi values latitude S2H OUT LAM lam raw Only for ambiguity method Output file name for computed lambda values longitude S2H NPOINTS 200 Only for schwabisch method the number of locations to compute the reference phase at different altitudes S2H DEGREE1D 2 Only for schwabisch method the degree of the 1d polynomial to fit reference phase through at every location S2H_NHEIGHTS s2h_degree1d 1 Only for schwabisch method the number of heights to evaluate the reference phase minimum is default S2H_DEGREE2D 5 108 Only for schwabisch method the degree of the 2d polynomial to fit 1d coefficients as function of location Example input c c comment SLANT 2 HEIGHT CONVERSION E S2H_METHOD schwabisch S2H NPOINTS 500 S2H_DEGREE1D 2 S2H_NHEIGHTS 3 S2H_DEGREE2D 5 SALMO DIET Outdata hei schw c S2H_METHOD ambiguity a SL OUT KINI Outdata hei ambi c S2H OUT LAM Outdata lam ambi e SL OUT IHE Outdata phi ambi e c S2H METHOD rodriguez e SiL OUT leur Outdata hei rodr
187. slave image M and S 2 If there is a reference phase R real values evaluate it at the master grid buffer and then subtract it in a complex way from slave S and store it in S matlab like notation pointwise multiplication notated by S S x cos R isin R 23 4 3 Compute complex interferogram and store it in M M M x conj S M conj S conj R 23 5 4 Multilook complex interferogram if requested Do not divide scale multilooked interferogram 5 Write this buffer to disk and start next one See Figure 23 2 for an explanation of the use of buffers in the implementation In this example the number of lines equals 17 The number of pixels equals P Suppose multilook factor in line direction mL 3 and in pixel direction mP 3 Furthermore after computing the available memory the number of lines of one buffer BL is maximum 7 BL is set to a multiple of mL BL 6 Now the number of fully filled buffers NB 17 6 2 The number of lines left in the last buffer equals 17We can compute something in this last buffer only for the first three lines so the last buffer is 3 lines long The first buffer is read master and slave image and computations are performed These results are written to disk Then next buffer etc If buffers 3 then resize the matrices for the computations The number of lines of the total result are L mL and P mP floored 76 Chapter 24 COMPREFPHA In this chapter the processing o
188. st oversampling the images itself and directly computing the correlation for a small number of shifts assuming initial offsets are known within a few pixels as we suspect that there may be an error introduced due to aliasing with the method that is implemented This method will be named oversample The actual computation of the transformation model 2d polynomial is done by the step COREGPM compu tation of coregistration parameters See also Samson 1996 17 1 Input Cards FC_METHOD magfft magspace oversample Select method for the computation Compute cross correlation based on magnitude im ages either in the space or the spectral domain Magnitude patches are zero meaned Method magfft is fast but patch size varies depending on shift Magspace keeps con stant patch size and shifts it over the master but is slower Computations are done in space domain Method oversample is best theoretically avoids aliasing of spectrum when magnitude is computed using FFTs FC_NWIN 400 The number of windows to be distributed over the total image if points are read from file FC IN POS then this card is ignored FC IN POS file name A ascii file with integer line pixel pairs coordinates in the original master coordinate system with locations where the windows should be placed After the last coordinate there should NOT be a EOL enter though Doris should ignore this The coordinates should be within the current overlap of
189. t gray steps are optional Obviously the modular structure of Doris enables the user to specify its own processing chain The most important intermediate products are indicated on the right hand side 1 2 General considerations and conventions In this section some important definitions are described that are used in the Doris software This will clarify the terminlogy used The general set up of the input file and output files is described as well After compilation with the Makefile the executable is named doris In this document therefor this name is used to refer to the executable The command line options are e doris ver return version number e doris h return help system call to shell script named helpdoris e doris file run use input in file default inputoptionsfile It is advised also to compile an executable doris debug with the Makefile This version is somewhat slower and more verbose but it can be used if something seems to go wrong with the normal executable and it is not clear what See also Annex B We advice to use the utility scripts to generate input files and to run the processor Conventions are e We use the term ines to refer to the azimuth direction slow time and pixels for the range direction fast time A pixel might also refer to an element which will be clear from the context In the source code we frequently use the term azimuth buffers and range blocks e Our convention
190. ted 3 3 1 Changes for X86 platforms Version 2 4 onward can be easily compiled on Linux systems For little endian machines like intel PC s this means the byte order is different Since the record length in the leader volume and data file are stored as B4 Abytes unsigned integers we had to use the function ntohl see the manual pages routines readvolume readleader readdat Also the SLC data itself in the datafile is stored as 2x 2B short signed integer byte data real part imaginary part real imaginary real imag We transform the data that is read from the data file with the function htons if X86PROCESSOR is defined see the Makefile or source code for more information 16 Chapter 4 M_PORBITS In this chapter the step M PORBITS is described This step should be run after READFILES because in that step the azimuth time is written to the master result file from from the SLC leader file We use the DEOS fortran program getorb for obtaining the precise orbits This program has to be installed separately see http www deos tudelft nl ers precorbs or Scharroo and Visser 1998 This step actually is only a system call to getorb and converts the output to a 4 column table secofday x y z It requires the Orbital Data Records ODR files to be in an archive directory The arclist file which should be downloaded together with the ODR archives has to be present in this directory The ODR files have to be untarred
191. ted also see the Makefile and matlib c We previously ran into some problems with the gnu g compiler version 2 95 2 and were not able to compile Doris explicit instantiation mechanism parsing of friend template functions bug in this version We solved this problem by putting all definitions of the matrix class in one file and compiling without first creating a matrix library Doris v2 4 and higher should therefor be more compiler platform independent For developers this is not quite so comfortable because the code is less transparent this way If compilation as described above gives problems try the following 1 Try make n processor o only echos command to the screen and paste this from the prompt to give you more direct control 2 Do not use Veclib and Lapack libraries even if you have them To exclude them define DEF4 and DEF5 in the Makefile uncomment the defines Now you will be making use of the internal routines which are based on numerical recipes routines 3 Set verbose flags for compilation v likely add to CFLAGS 4 Try to compile object code for a individual source files e g make processor o and if successful try others make swobjs does all 5 f all o files are compiled correctly use make or make doris to link them 6 Try another compiler first make clean B 5 3 Some notes on installation on SGI This message was posted to the doris users email list It may be of help if you are installing D
192. tern Return help for search pattern Calls the shell script helpdoris doris q Return random quote This option can be used to add random quotes to your mail Make an alias for e g elm or pine mailprograms alias elm doris q signature elm Then the next time elm is called it first creates a signature file in your home directory which is appended to your mail message doris c Return copyright notice doris inputfile Run doris with the input specified in inputfile defaults to inputoptionsfile For convenience we have developed a simple csh script named run in the BIN directory that can generate a template input for you and that works as a shell for the actual processing You can run the Doris software by for example assume you know the vi editor 1 Make a directory mkdir data kampes Testdoris Go to that directory cd data kampes Testdoris 2 Copy the run file for editing cp nome Doris BIN run Edit the run file for bin directories and remarks vi run Take a look at the help run h Generate input templates run g 3 Edit the generated input file run e1 Run the first step run s1 127 4 View output stout run v1 View output result file run r1 View output logfile run r4 5 Process next step run e2 run s2 run v2 etc B 4 Viewing the results of Doris Annex C describes a number of tools to visualize the output of Doris The parameter files log file etc can
193. th fixed reference system 1 1 1 Processing order The processing order is not restricted but in general the output of a step is the input for the next The following flowchart Figure 1 2 shows the general processing order The black processing steps are obligatory to obtain a geocoded interferogram as end product The dark gray steps are highly recommended whereas the light gray steps are optional Obviously the modular structure of Doris enables the user to specify its own processing chain The most important intermediate products are indicated on the right hand side Processing steps Products 3 M_READFILES 9S READFILES 4M PORBITS 10S PORBITS 5M CROP 11S CROP Crops 6 M_SIMAMP 7 M_TIMING Oversampled crops 13 COARSEORB 14 COARSECORR 15 M_FILTAZI E 16 S FILTAZI 17 FINE 18 RELTIMING 19 DEMASSIST 20 COREGPM 21 RESAMPLE Resampled slave 22 FILTRANGE 23 INTERFERO Interferogram 24 COMPREFPHA Interferogram 25 SUBTRREFPHA subtracted reference phase 26 COMPREFDEM Interferogram 27 SUBTRREFDEM subtracted reference DEM 28 COHERENCE 3 Coherence image 29 FILTPHASE S0 UNWHAP Unwrapped interferogram EE Obligatory step al Recommended step 32 SLANT2H SS GEOCODE Geocoded interferogram Figure 1 2 Processing flowchart for Doris The black processing steps are obligatory to obtain a geocoded interferogram as end product The dark gray steps are highly recommended whereas the ligh
194. the satellite height yu the location angle of the co registered resolution cell in the interferogram angle between two vectors x and y by Z rz y B d M S D 1 B 4 M P d S P D 2 Now the perpendicular baseline has to be computed The definition states that B is positive if the slave satellite is to the right of the slant range line of the master Which yields for a mountain an increasing phase from foot to summit We had some trouble finding a simple expression to find out the correct sign but at the moment we do something like the following B B B D 3 rn M P D 4 n S P D 5 y Z P F D 6 2 Z P Fa D 7 _ 1 jm lt y2 sign i n wl D 8 145 Table D 1 Conversion between baseline representations note that the four quadrant arctangent should be used Bn By B a B1 Bj Ba Ba B Bcosa Bj B cos0 Bj sind Il B Bsina B B sin Bj cos0 B a a arctan B Bi a arctan B B_ B yB B B wy By B1 B B B B sin0 B cos By Bsin 0 a Il Il Il B B c0s0 B sin B Bcos 0 0o 9 Z M T1 a 0 arctan 2 B B1 D Bj B cos a D 1 B Bsin a D D 3 Interferogram The phase for a certain pixel in a single SLC image is defined as fap D 14 The complex interferogram minus reference phase is defined as I M S R D 15 Where denotes the complex conjugated a 2 2
195. thod uses the geometry to compute sin 0 a There are two errors in it for now First H is not computed exact Second the baseline parameters are not computed exact per line P is evaluated at reference surface and S according to that position This means that if the orbits are not parallel the point S is not computed correctly which introduces errors in the baseline computation The co registration model is better used for that Known r1 range to M P position M B baseline baseline orientation with regard to equals our def of alpha Compute theta angle state look with formulas exact no iterations H height sat above some surface how to compute this exact The following we derived for our baseline definition beta angle 2 1 P 1 counterclockwise B Z 2 1 P 1 32 11 cosB cos 7 a 8 32 12 sin 0 a 32 13 r By r B 204 B cos B 82 14 111 r By r B sin 0 a 3B5 32 15 Bee 32 16 So theta can be solved for exact with these formulas Note B 2_p2_ B 0 a arcsin eens 32 17 0 arcsing a V 0 7 arcsinr a 32 18 did not find an efficient way to always use the right expression yet Now use the fact that theta is about 20 degrees but it should be possible to find out quadrants Compute H from known theta position Master rho1 r1 In triangle 1 P 0 three Start with cosine law for line across theta p p p r 2rhoi
196. tion NaN not a number correlation values which are mainly due to no data regions at the border of the scene are eliminated automatically but the whole list of windows are kept in the doris log file Coherence output is properly scaled between 0 to 1 range with new option r of cpxfiddle see Chapter C 2 3 123 Annex B Installation In this annex the installation of Doris is described To properly compile the Doris software you might have to edit the Makefile Set the compiler CC compiler flags CFLAGS library path LFLAGS and defines comment DEF4 and DEF5 for VECLIB LAPACK library usage Use DEF7 _X86PROCESSOR_ if you have a little endian machine We have written a simple script configure to help generate a user defined Makefile which is present in the Doris distribution as well as a template Makefile that can be edited to your likings if the script fails We recommend compiling 2 versions of the Doris software An optimal version for operational processing and a more verbose debug version that is only to be used if the optimal version exits with an unexplained error Then repeat processing with debug version track down routine etc Compilation of these two versions is best done by running make two times first with CFLAGS CFLAG SOPT then after make clean with CFLAGS CFLAGSDEBUG This is clearly described in the Makefile and in the Makefile generator We have successfully compiled Doris with GNU g 2 x
197. tion is written 18 3 Implementation The observation equations are given by the zero degree polynomial model y A x y 1 Y2 1 ai 0 p 0 18 1 UN 1 Where y contains the observed offsets in a certain direction Qip denotes the unknown coefficients of the polynomial The least squares parameter solution is given by ATQ y ATQ Az Na 18 2 Where Oo is the diagonal covariance matrix of the observations The coefficients are estimated by factorization of the matrix N The inverse of matrix N is also computed to check the solution stability and to compute some statistics A check number is given max abs N N I that gives a hint on the stability of the solution 49 Chapter 19 DEMASSIST In this chapter the processing of step DEMASSIST is described Here the slave image is coregistered to the master image based on a DEM For each pixel of the master image the corresponding real valued coordinate in the slave image is computed For this step a DEM is required e g obtained by the SRTM mission The Doris distribution contains the utility construct_dem sh do download and prepare SRTM data see Section C 2 12 The utility also outputs a figure of the final result ps format and the lines needed for the Doris input file The DEM assisted coregistration is not dependent on the correlation between the master and the slave image Coregistration errors due to bad distribution or lack of useful corr
198. tive should be performed after the resampling of the slave to the master grid because the fringe frequency is estimated from the interferogram that is temporary computed It is performed simulataneous for the master and slave image 22 1 Input Cards RF METHOD adaptive porbits Method selector for range filtering Either adaptive recommended or based on the precise orbits RF FFTLENGTH 64 For method porbits and adaptive For method porbits length of block in range direction 512 or 1024 default for this method advised For method adaptive Length of window for adaptive method A peak is estimated for parts of this length RF OVERLAP 0 For method adaptive Overlap between input buffers in range direction RF HAMMING 0 75 For method porbits and adaptive Weight for hamming filter 1 is rect 66 RF_SLOPE 0 For method porbits Terrain slope in degrees Positive slope is towards radar A slope equal to the viewing angle implies total filtering RF NLMEAN 15 For method adaptive Take walking mean over RF_NLMEAN lines to reduce noise for peak estimation Has to be odd Compare with periodogram RF_THRESHOLD 5 For method adaptive Threshold on SNR of peak estimation to perform range filtering RF_OVERSAMPLE 2 For method adaptive Oversample master and slave with this factor before computing the complex interferogram for peak estimation This factor has to be a power of 2 2 is default to be able to estimate the peak for freque
199. to account for most effects in normal images The 2d polynomial has the form 53 Location_windows wtests Range 399 798 1197 1596 1995 2394 2793 3192 3591 3990 14990 Hee Hee 1e Hot 1 I L I 14990 Porti 478 4480 a81 4 483 467 469 472 3 M54 M55 H56 H57 M58 459 13491 M42 48 4444 445 M46 0947 448 L13491 34 436 420 421 441422 4423 425 4426 p407 4 409 410 aen 113 B98 4401 4403 11992 4 884 4885 4 886 4 B87 B88 B89 4 890 11992 4 874 4 B76 877 B78 4 79 67 lt 869 4857 10493 Seas B46 10493 p33 8994 4 et in a T 8994 iz 4 B74 2 hp65 4 252 253 L E 94 REI aes 745 E 2 831 2 4 18 05 09 5996 TR ri ener 5996 70 4 471 74 di 142 4497 46 150 gt 4497 att 37 438 39 4423 4427 11 42 4413 Sus 02 h 2998 A o hi b 2998 4 80 4 3 p6 1499 ua F 1499 81 0 ez 43 fe Se 4 T T T T T T T T T 399 798 1197 1 sog 1995 2394 27 3192 3591 3990 Ee 16 0556200 Range Figure 20 1 Plot produced by plotcpm for the first run The estimated offsets are plotted here normalized together with a 90 degrees rotated w test as ellipses Azimuth direction Correlation 0 0 0 1 0 2 0 3 0 4 0 5 0 6 0 7 0 8 0 9 abs e abs e 00 01 02 03 04 05 06 07 08 09 E Correlation Figure 20 2 Plot produced by plotcpm for the first run The absolute error estimated offsets minus observed offsets are plotted fo
200. tory LOGFILE log out Output filename for the logfile M RESFILE master result out Output filename for the master result file S RESFILE slave result out Output filename for the slave result file RESFILE interferogram out Output filename for the products result file ORB INTERP POLYFIT DEGREE j SPLINE 11 Orbit interpolation method Defaults to a polynomial of degree numberofdatapoints 1 but smaller than degree 5 order 5 Optionally the DEGREE can be given j numberofdatapoints 1 The x y z are independently interpolated the velocities are estimated from the position If method SPLINE is selected natural cubic splines are used This may be inaccurate if there are only a few orbit datapoints default interpola tion not approximation for polyfit is used to go smoothly through the datapoints since the points do probably not contain noise since they are already the result of an orbit propa gator somewhere It is not adviced really to use a DEGREE smaller than the maximum possible except if it gets too large to avoid oscillations DUMPBASELINE 00 Dump the baseline parameters for a grid of O lines by 0 pixels as INFO to stdout The baseline is only evaluated after the orbits are known The perpendicular baseline to the reference ellipsoid is also computed as a 2D polynomial of degree 1 And also theta as function of azimuth line and range though it hardly varies over azimuth HEIGHT 0 0 Average terrain height ab
201. tput window You can make a cutout of the image with this card Hi values larger than the size of the image are reset to the maximum If card omitted it defaults to total image line pixel 1 refers to the the first line pixel M DBOW GEO lat_0 lon_0 height width Master output window Alternative to and overrides normal DBOW card You can make a cutout of the image with this card latitude of the center pixel of the desired crop longitude in decimal degrees WGS84 system of orbit then height width in pixels For approximately square areas heights should be a factor 5 of width for ERS Example input cards for this step e comment CROP 20 a M_CROP_IN cdrom SCENE1 DAT_01 001 M_CROP_OUT Output 21066 raw c M_DBOW iL SOO i 10010 linelow hi pixellow hi M_DBOW_GEO 52 6734 5 4342 10000 2000 lat_O deg lon 0 height width pix 5 2 Output Description The process control flag at the start of the result file is switched to 1 at successful exit CIEOIO S il The output section in the result file will resemble the following kk kk kk kk Sk Sk Sk Sk KK KK KK KK KK KK KK KK KK KKK KKK KKK ARK KKK kk kk Ck Ck Sk Ck Pk A A A k kx ko ARA 5 SIEGE eso master HA RARA Fee KE ee ee oe AA AAA oe to ie AR AAA toe Ak at ake AAN A AR AAA te eae Xo ee ak ae kc t ox Data_output_file Output 21066 raw Data_output_format complex_short First line w r t original master 1 Last line w r t original master 5000
202. ts c void dinsar 138 writearg readinput c void writearg const Type argument 139 readinput readinput c void readinput 140 checkgeneral readinput c void checkgeneral 141 checkreadfiles readinput c void checkreadfiles 142 checkcrop readinput c void checkcrop 143 checkporbits readinput c void checkporbits 144 checkslant2h readinput c void checkslant2h 145 checkunwrap readinput c void checkunwrap 146 checkgeocode readinput c void checkgeocode 147 checkcoarsecorr readinput c void checkcoarsecorr 148 checkfine readinput c void checkfine 149 checkcoregpm readinput c void checkcoregpm 150 checkcomprefpha readinput c void checkcomprefpha 151 checksubtrrefpha readinput c void checksubtrrefpha 152 checkresample readinput c void checkresample 153 checkinterfero readinput c void checkinterfero 154 checkcoherence readinput c void checkcoherence 136 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 checkcomprefdem checksubtrrefdem checkfiltrange checkdinsar checkfiltphase checkfiltazi setunspecified specified flatearth radarcodedem fillslcimage updateslcimage readinput readinput readinput readinput readinput readinput readinput readinput h referencephase referencephase
203. turned on ELE Ail 2 I In the master result file a section is added with the new file name okokok oko ck ko ko AEREA A ARA ERE RARA ARE kokok kok kok k ck kok ko kk kok AA ARA EKER KEK kk k a Sheet seal lle ean B KEKKKKKKKKKKKKKKKKKKKKKKKKKK KK KKK KKK KKK ckckckckckckckckckckckckck ck ck ck ck ck ck ckckck ck ck kk kk Input file Outdata 1393 raw Data output file Outdata 1393 azifilt Data output format complex real4 Elis to aga ai X SE S 3 Ti Last line w r t original master 3500 First pixel w r t original master Il Last pixel w r t original master 500 XC se eK X Fe He a IA Ae ee Fe He a e AA A ARA to he ct oe c He a ak RAR oe CS koe Xo eRe Xo ok ak ake ko ox End filt azi NORMAL Ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck ck Ck ck ck ck Ck ck ck ck Ck ck ck Ck Ck Ck Ck Ck CK CK CK CC CK CIS kk Sk Sk Sk I ke S ke kk Sk kc k ck ck ckckok A file mph is created for the master with filtered spectrum Figure 15 1 demonstrates the filter for 2 images 15 3 Implementation For each buffer of AF BLOCKSIZE lines and width pixels do taking care of AF OVERLAP e Take 1DFFT in azimuth direction over the columns 40 MASTER SLAVE ER T T T 150007 15000r 10000rF 10000F 5000 5000F pla all MMC UR i A M ol TIENES LL H 1 T f 800 600 400 200 0 200 400 600 800 800 600 400 200 0 200 400 600 800 filter for master red is compo
204. ures 31 1 31 2 and 31 3 give example of output 31 3 Implementation See also Zebker et al 1994 See figure 31 4 Simple equations for topo pair no deformation no atmosphere no other errors r r2 0 00 00 31 1 By di T 31 2 B Bcos 0 a B cos a 0 31 3 By Bsin 0 a Bsin a 0 31 4 The baseline components for points on the reference ellipsoid 0 are B1 9 Bio B cos 0o a 31 5 Bi 0 Bio Bsin o a 31 6 The true phase of the interferogram is 4 jc XB 81 7 103 dou MIT Be bid eec mal dm nm Grn Figure 31 1 Phase of complex topography interferogram flat earth corrected cropped The area is dead sea Israel Temporal baseline is 1 day tandem The perpendicular baseline is approximately 105 meters This interferogram has been coregistered on the defo pair 31 2 by tricking Doris TH H angr r eae DE Figure 31 2 Phase of complex deformation interferogram flat earth corrected cropped The area is dead sea Israel Temporal baseline is 28 months day The perpendicular baseline is approximately 30 meters 104 a Tem diem ill Tip oor i Figure 31 3 Phase of differential complex interferogram result of step DINSAR cropped The area is dead sea Israel The topography is removed from the original interferogram Figure 31 2 by scaling the topography interferogram Figure 31 1 The perpendicular b
205. ut Cards 30 2 Output Description ouo ete hoe eee een pug Y Y y re RSS 31 DINSAR 31 1 Input Cards 31 2 Output Description oko bebe EP ER o oon DoD mh om Rom Y ELA EUR ow E D Ses 31 3 Implementation 31 3 1 Algorithm 32 SLANT2H 32 1 Input Cards 32 2 Output Description ns 32 3 Implementation 32 3 1 Method ambiguity s r sas e sow aaa a ano a aa a a A a a a a oo i a 32 3 2 Method rodriguez es 32 3 3 Method schwabisch es 32 4 Comparison ofthe methods s s s a sos aaraa a llle 33 GEOCODE 33 1 Input Cards sn ea mae Se A ida a eee EUR Ree BS 33 2 Output Description ee A 33 3 Implementation ee ee A What s new A 1 A 2 Version Pu ee es oo Race ae Eh aha c ds ie Version AOT xx a a a a aa eee os xe hk REE a aaa S E a anae RR OR RAS EE em we AY B Installation B 1 B 2 B 3 B 4 B 5 B 6 B 7 Installation Of Doris 4 ge beh eee EUER a RC Gh Ge RUEME Bc ee B 1 1 Installation of the Doriscore 0 o lll B 1 2 Installation of the SARtools 0 a B 1 3 Installation of the ENVISAT tools 2l B 1 4 Installation of the TERRASAR X reader lll ln B 1 5 Starting DOSa es he e x mam mom Pe A ee eee ere UR HR IR IQ Re a ARN B 1 6 Installation of utility scripts 2 2 2 2 e Additional programs i a as e s ao hm oe Running the Doris software ee Viewing the results of Doris ee
206. will write this in the Product Type Specifier field RSAT must be tested problems may be orbit data In later steps the Product field is read and the CROP step uses the appropriate function automatically Envisat ERS JERS RSAT ATLANTIS M IN DAT filename The filename of the SLC data file This is the only file required for method ASAR EN VISAT M IN LEA filename The filename of the SLC leader file Not used for method ASAR ENVISAT M IN VOL filename The filename of the SLC volume file Not used for method ASAR ENVISAT and TSX TERRASAR X M IN NULL filename The filename of the SLC null file This may be a dummy name since it is not used 14 An example of the input cards for this step is given below This example can be inserted in the general cards described in Section 2 2 E comment READFILES e M IN METHOD ERS M IN VOL cdrom scenel vdf dat 001 name of volumefile M IN LEA cdrom scene1 lea 01 001 name of leaderfil M IN NULL dummy name of nullfile M IN DAT cdrom scenel dat 01 001 name of datafile 3 2 Output Description The process control flag at the start of the master result file is switched to 1 at successful exit readfiles il Example of output of this step in master result file The positioning data of the platform from the leader file has been deleted in this example This happens automatically after step M PORBITS getting the precise orbits Thi

Download Pdf Manuals

image

Related Search

Related Contents

Yamaha RX-V1300 6.1 Channels Receiver  User Manual FavorPrep Soil DNA Isolation Midi Kit  コンビニ収納代行ソフト「コンペイ君」取扱説明書  FPGA Interface Board User Manual    Origin Storage i4500R 1U 4 Bay Rackmount iSCSI SAN 1000GB  Neorest® II Lavatory Neorest® II Lavabo Neorest® II Lavabo  Cargar  Category 10/セキュリティ(P.185~188)  ActiveJet AB-M1100R  

Copyright © All rights reserved.
Failed to retrieve file