Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

GACOS: Error due to lat/lon and UTM coordinates for HYP3 products #664

Open
Sugeradi985 opened this issue Sep 18, 2021 · 3 comments
Open

Comments

@Sugeradi985
Copy link

Description of the problem
An error occurred when I used Mintpy to do atmospheric correction (GACOS) on the interferometric pairs derived from HYP3 sdk.

I think it is the coordinate system differences that mattered but not sure, and GACOS data (WGS84) cannot be read normally when processing interferometric pairs (WGS84 UTM50N) derived from HYP3 sdk.

******************** step - correct_troposphere ********************
Input data seems to be geocoded. Lookup file not needed.
tropospheric delay correction with gacos approach

tropo_gacos.py -f /mnt/hyp3_mintpy_asf/bj_142_126/newds/mintpy/timeseries.h5 -g /mnt/hyp3_mintpy_asf/bj_142_126/newds/mintpy/inputs/geometryGeo.h5 -o /mnt/hyp3_mintpy_asf/bj_142_126/newds/mintpy/timeseries_GACOS.h5 --dir ./GACOS
Use GACOS products at directory: /mnt/hyp3_mintpy_asf/bj_142_126/newds/mintpy/GACOS
update mode: ON
output file: /mnt/hyp3_mintpy_asf/bj_142_126/newds/mintpy/inputs/GACOS.h5
1) output file either do NOT exist or is NOT newer than all ZTD files.
run or skip: run
--------------------------------------------------
create HDF5 file: /mnt/hyp3_mintpy_asf/bj_142_126/newds/mintpy/inputs/GACOS.h5 with w mode
create dataset  : date       of |S8                       in size of (120,)               with compression = None
create dataset  : timeseries of <class 'numpy.float32'>   in size of (120, 2799, 3591)    with compression = None
close  HDF5 file: /mnt/hyp3_mintpy_asf/bj_142_126/newds/mintpy/inputs/GACOS.h5
read incidenceAngle from file: /mnt/hyp3_mintpy_asf/bj_142_126/newds/mintpy/inputs/geometryGeo.h5

Full error message

Traceback (most recent call last):
  File "/home/melon/tools/MintPy/mintpy/smallbaselineApp.py", line 1279, in <module>
    main(sys.argv[1:])
  File "/home/melon/tools/MintPy/mintpy/smallbaselineApp.py", line 1261, in main
    app.run(steps=inps.runSteps)
  File "/home/melon/tools/MintPy/mintpy/smallbaselineApp.py", line 1054, in run
    self.run_tropospheric_delay_correction(sname)
  File "/home/melon/tools/MintPy/mintpy/smallbaselineApp.py", line 791, in run_tropospheric_delay_correction
    mintpy.tropo_gacos.main(iargs)
  File "/home/melon/tools/MintPy/mintpy/tropo_gacos.py", line 361, in main
    calculate_delay_timeseries(tropo_file=inps.tropo_file,
  File "/home/melon/tools/MintPy/mintpy/tropo_gacos.py", line 306, in calculate_delay_timeseries
    delay = get_delay_geo(ztd_file, atr, cos_inc_angle)
  File "/home/melon/tools/MintPy/mintpy/tropo_gacos.py", line 129, in get_delay_geo
    delay = resize(delay, (length, width),
  File "/home/melon/tools/miniconda3/envs/insar/lib/python3.8/site-packages/skimage/transform/_warps.py", line 176, in resize
    out = warp(image, tform, output_shape=output_shape, order=order,
  File "/home/melon/tools/miniconda3/envs/insar/lib/python3.8/site-packages/skimage/transform/_warps.py", line 827, in warp
    raise ValueError("Cannot warp empty image with dimensions",
ValueError: ('Cannot warp empty image with dimensions', (0, 0))

System information

  • Operating system: Linux-Ubuntu 18.04.5
  • Python environment: conda
  • Version of MintPy: MintPy release version v1.3.1
@welcome
Copy link

welcome bot commented Sep 18, 2021

👋 Thanks for opening your first issue here! Please make sure you filled out the template with as much detail as possible. We appreciate that you took the time to contribute!
Make sure you read our contributing guidelines

@yunjunz yunjunz changed the title Error occured during atmospheric correction (based on GACOS) of HYP3-derived interferometric pairs GACOS: Error due to lat/lon and UTM coordinates for HYP3 products Sep 18, 2021
@yunjunz yunjunz added this to the UTM coordinate support milestone Sep 18, 2021
@yunjunz
Copy link
Member

yunjunz commented Sep 18, 2021

Thank you @Sugeradi985 for reporting this bug. I agree with you that it's due to the different coordinate systems between GACOS and HyP3 products. I would suggest the following two solutions:

  1. crop the GACOS products based on the SNWE of HyP3; then resize the 2D matrix. The step by step procedure would be:
  • grab the four corner coordinates of the HyP3 products in UTM;
  • convert them to lat/lon using ut.to_latlon();
  • use these lat/lon to read/crop the GACOS products;
  • resize the cropped GACOS products to the same size as the HyP3 products using skimage.transform.resize(). Check similar usage in tropo_gacos.py#L128 for reference.
  1. resample the 2D matrix from lat/lon coordinates to UTM coordinates.
  • calculate the pixel-wised coordinates of HyP3 products in lat/lon using ut.to_latlon()
  • resample the GACOS 2D matrix from its original lat/lon grid to the lat/lon coordinates calculated above.

Option 1 is used in tropo_pyaps3.py currently. Option 2 is definitely the right way, but can be slower and slightly more complex than option 1. I am not 100% clear about the difference between the two.

@mosyhey
Copy link

mosyhey commented Jun 18, 2022

You need to add this piece of code to tropo_gacos.py, line 120 (after geo_box part):

    ########## New lines for converting utm to lat/lon (https://github.com/insarlab/MintPy/issues/664)
    strProcessor = atr['PROCESSOR']
    if strProcessor == 'hyp3':
        import utm #Install utm as "pip install utm" if not installed
        strUtmZone = atr['UTM_ZONE'] #Read from metadate, for example: "11N"
        intUtmZone = int(strUtmZone[0:-1]) #In this example: 11
        strHemisphere = strUtmZone[-1] #In this example: "N"
        north, west = utm.to_latlon(geo_box[0], geo_box[1], intUtmZone, strHemisphere)
        south, east = utm.to_latlon(geo_box[2], geo_box[3], intUtmZone, strHemisphere)
        geo_box = (west, north, east, south)
    ##########

Be aware that there is a very small and negligible rotation in converting a shape from UTM into Lat/Lon coordinate system.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Development

No branches or pull requests

3 participants