insar_vs_gps function problem

[Melanie Ayala]

I'm attempting to use the insar_vs_gps function in MintPy. I am using it to evaluate the accuracy of vertical displacement estimates generated from Sentinel 1 on hyp3. This processing chain generates the following files in the inputs folder: ERA.h5, geometryGeo.h5, ifgramStack.h5, mintpy_config.txt, smallbaselineApp.cfg. The smallbaselineApp.cfg file is attached. 

I am attempting to follow this notebook to understand usage of the insar_vs_gps function, and I am hitting a wall with my understanding of how to manipulate mintpy geometry files. The original code:

# GPS site_names = ['GV{:02d}'.format(i) for i in [1] + np.arange(3,11).tolist()] gps_dir = os.path.expanduser('~/data/Galapagos/GPS') ref_site = 'GV01' gps_loc_file = os.path.join(gps_dir, 'loc.txt') start_date = '20141101' end_date = '20180625' # InSAR proj_dir = os.path.expanduser('~/insarlab/Galapagos/GalapagosSenDT128/mintpy') out_file = 'mintpy.pickle' ts_file       = os.path.join(proj_dir, 'timeseries_ECMWF_ramp_demErr.h5') geom_file     = os.path.join(proj_dir, 'inputs/geometryRadar.h5') temp_coh_file = os.path.join(proj_dir, 'temporalCoherence.h5') if not os.path.isfile(out_file):     obj = insar_vs_gps(ts_file, geom_file, temp_coh_file,                        site_names, gps_dir, ref_site=ref_site,                        start_date=start_date, end_date=end_date)     obj.open()     with open(out_file, 'wb') as f:         pickle.dump(obj.ds, f)     print('saved to pickle file {}.'.format(out_file))     ds = dict(obj.ds) else:     with open(out_file, 'rb') as f:         ds = pickle.load(f)     print('read from pickle file {}.'.format(out_file))

Since my dataset is in geographic coordinates instead of radar coordinates (has a geometryGeo.h5 file instead of a geometryRadar.h5 file), I attempted to use the geometryGeo.h5 file instead of geometryRadar.h5 file as the geom_file input to the function:

# InSAR proj_dir = os.getcwd() out_file = 'mintpy.pickle' ts_file       = os.path.join(proj_dir, 'timeseries.h5') geom_file     = os.path.join(proj_dir, 'inputs/geometryGeo.h5') temp_coh_file = os.path.join(proj_dir, 'temporalCoherence.h5')

and that results in the following error message

--------------------------------------------------------------------------- ValueError Traceback (most recent call last) Cell In[28], line 3 1 if not os.path.isfile(out_file): 2 obj = insar_vs_gps(ts_file, geom_file, temp_coh_file, site_names,start_date=start_date, end_date=end_date) ----> 3 obj.open() 4 with open(out_file, 'wb') as f: 5 pickle.dump(obj.ds, f) File /opt/homebrew/anaconda3/envs/coastal-py/lib/python3.10/site-packages/mintpy/objects/insar_vs_gps.py:102, in insar_vs_gps.open(self) 99 raise ValueError(msg) 101 self.read_gps() --> 102 self.read_insar() 103 self.calculate_rmse() 104 return File /opt/homebrew/anaconda3/envs/coastal-py/lib/python3.10/site-packages/mintpy/objects/insar_vs_gps.py:141, in insar_vs_gps.read_insar(self) 139 geo_box = (min(lons), max(lats), max(lons), min(lats)) #(W, N, E, S) 140 pix_box = coord.bbox_geo2radar(geo_box) #(400, 1450, 550, 1600) --> 141 src_lat = readfile.read(self.geom_file, datasetName='latitude', box=pix_box)[0].reshape(-1,1) 142 src_lon = readfile.read(self.geom_file, datasetName='longitude', box=pix_box)[0].reshape(-1,1) 143 src_pts = np.hstack((src_lat, src_lon)) File /opt/homebrew/anaconda3/envs/coastal-py/lib/python3.10/site-packages/mintpy/utils/readfile.py:296, in read(fname, box, datasetName, print_msg, xstep, ystep, data_type, no_data_values) 294 fext = os.path.splitext(os.path.basename(fname))[1].lower()

...

--> 375 raise ValueError(f'input dataset {datasetName} not found in file {fname}') 377 # output size for >=2D dataset if x/ystep > 1 378 xsize = int((box[2] - box[0]) / xstep) ValueError: input dataset ['latitude'] not found in file /Users/melanieayala/Documents/Code/inputs/geometryGeo.h5

OS: MacOs 13.3 M2 MacBook Air

Mintpy installation method: conda

Mintpy version: MintPy version .0.26, date 2023-06-27

My three questions:

1. Is there a mintpy function that converts a geometeryGeo.h5 file to a geometryRadar.h5 file, and vise versa? If so, what is the function, and can I get a usage example? 
2. Can I run the insar_vs_gps function with the geometeryGeo.h5? If so, can I get a usage example?
3. How would I modify the smallbaselineApp.cfg file to produce a geometryRadar.h5 file in addition to or instead of the geometeryGeo.h5 file?

[Zhang Yunjun]

You may find a relevant discussion here (https://github.com/insarlab/MintPy/discussions/1010). The insar_vs_gps module is not as generic as the rest of mintpy. There are two issues with using it for the hyp3+mintpy workflow: 1) the missing "latitude" dataset in the geometry file, which leads to the error you have; and 2) the UTM coordinates, which are not supported yet. Unfortunately, there are no ready-to-use solutions at the moment. Scott Staniewicz from JPL has issued a PR (https://github.com/insarlab/MintPy/pull/1079), which aims to fix the 2nd issue, but it's not ready to be merged yet.

Thus, I would suggest to NOT using the insar_vs_gps module for hyp3+mintpy at the moment, but instead, write your own code for this comparison purpose, following a similar procedure as insar_vs_gps module:

1. download and read the GNSS data (from UNR using the more generic gps module)

2. convert the GNSS lat/lon into UTM (using mintpy.utils.utils0.latlon2utm)

3. interpolate InSAR data at the GNSS location; OR select multiple InSAR pixels around the GNSS location (two flavors of comparisons that I am aware of)

4. find the common dates between InSAR and GNSS time series. 

I hope this is helpful,

Yunjun

[Zhang Yunjun]

Since this is a frequently asked question, I would like a make a maintained notebook example in the future. Could you open a issue on the mintpy-tutorial repo (https://github.com/insarlab/MintPy-tutorial) to remind us?

[DU SEN]

Hi Yunjun,

I am using gamma+mintpy, and I find that in the load_data.py function,  line 451, latitude/longitude is removed for gamma and roipac. I wander if this is related with the file inputs/geometryRadar.h5. Because I need latitude/longitude in  eometryRadar.h5 for my small program. Do you have any idea how can I add these two dataset to  inputs/geometryRadar.h5? 

Thanks