#!/usr/bin/env python
# coding: utf-8

# # Merge HyP3 ISCE2 burst InSAR products
# 
# This notebook demonstrates how to use the `merge_tops_burst` workflow of the HyP3-ISCE2 plugin. This workflow merges multiple burst InSAR Products and takes a directory that includes multiple HyP3-ISCE2 Burst InSAR Products as its input. These input products can be created by the HyP3-ISCE2 on-demand service. To learn how to create HyP3 burst InSAR Products, check out our [hyp3_isce2_burst_stack_for_ts_analysis](https://github.com/ASFHyP3/hyp3-docs/blob/main/docs/tutorials/hyp3_isce2_burst_stack_for_ts_analysis.ipynb) notebook.
# 
# 
# **Note:** This notebook does assume you have some familiarity with InSAR processing with HyP3 already, and is a minimal example without much context or explanations. If you're new to InSAR and HyP3, we suggest checking out the following resources. Note that some of these resources may be specific to our InSAR GAMMA products, so you may need to adapt them for use with our ISCE2-based burst InSAR products.
# 
# * Our [Burst Data Download Story Map](https://storymaps.arcgis.com/stories/88c8fe67933340779eddef212d76b8b8)
# 
# * Our [product guide](https://hyp3-docs.asf.alaska.edu/guides/burst_insar_product_guide/) for burst InSAR products
# 
# * Our [GitHub repository](https://github.com/asfhyp3/hyp3-isce2) containg the workflow used to create burst InSAR products
# 
# * Our [InSAR on Demand Story Map](https://storymaps.arcgis.com/stories/68a8a3253900411185ae9eb6bb5283d3)

# ## 0. Initial Setup
# To run this notebook, you will need a local copy of the HyP3-ISCE2 GitHub repository and to set up a conda environment with the required dependencies. In your terminal, you can do this with the following commands:
# 
# ```shell
# git clone https://github.com/ASFHyP3/hyp3-isce2.git
# cd hyp3-isce2
# mamba env create -f environment.yml
# mamba activate hyp3-isce2
# python -m pip install -e .
# mamba install -c conda-forge pandas jupyter ipympl
# jupyter notebook hyp3_isce2_burst_merge.ipynb
# ```

# This workflow requires an Earth Data Cloud login. If you haven't yet, you can make an [account for free](https://urs.earthdata.nasa.gov/users/new) and set up a `.netrc` file in your home directory with your personal username and password.
# 
# ```shell
# echo ‘machine urs.earthdata.nasa.gov login $USERNAME password $PASSWORD’ >> ~/.netrc 
# ```

# ## 1. Create burst InSAR prodcuts to merge
# 
# Before using the [HyP3-ISCE2 merge burst workflow](https://hyp3-docs.asf.alaska.edu/guides/burst_insar_product_guide/#merge-sentinel-1-burst-insar-products), we must create burst InSAR products that are merge-compatible. This means that they must:
# - Have the same reference and secondary dates
# - Have the same polarization
# - Have the same multilooking
# - Be from the same relative orbit
# - Be contiguous
# 
# In this section, we'll create two such burst InSAR products and download them to a local directory.

# In[ ]:


from pathlib import Path


reference_granule1 = 'S1_136231_IW2_20200604T022312_VV_7C85-BURST'
secondary_granule1 = 'S1_136231_IW2_20200616T022313_VV_5D11-BURST'

reference_granule2 = 'S1_136232_IW2_20200604T022315_VV_7C85-BURST'
secondary_granule2 = 'S1_136232_IW2_20200616T022316_VV_5D11-BURST'

project_name = 'merge_demo'
current_dir = Path.cwd() # Make sure we preserve this in case we want to navigate back
work_dir = Path.cwd() / project_name
data_dir = work_dir / 'data'
data_dir.mkdir(parents=True, exist_ok=True)


# In[ ]:


import hyp3_sdk as sdk

hyp3 = sdk.HyP3(prompt=True)

jobs = sdk.Batch()
for reference, secondary in [(reference_granule1, secondary_granule1), (reference_granule2, secondary_granule2)]:
    jobs += hyp3.submit_insar_isce_burst_job(
        granule1 = reference,
        granule2 = secondary, 
        apply_water_mask = False,
        name = project_name,
        looks = '20x4'
    )


# In[ ]:


jobs = hyp3.watch(jobs)


# In[ ]:


from datetime import datetime

now = datetime.now()
start_of_today = datetime(now.year, now.month, now.day)

jobs = hyp3.find_jobs(name=project_name, start=start_of_today)


# In[ ]:


insar_products = jobs.download_files(data_dir)
insar_products = [sdk.util.extract_zipped_product(ii) for ii in insar_products]


# ## 2. Merge the products using hyp3-isce2 merge_tops_burst workflow
# 
# Now we have all the data we need and can merge these two burst InSAR products! You would typically run the command below on the command line, but we'll run it through the Jupyter Notebook here.

# In[ ]:


cd $work_dir


# In[ ]:


# Use this command to actually view the options for the merge_tops_burst workflow
get_ipython().system('python -m hyp3_isce2 ++process merge_tops_bursts --help')


# In[ ]:


# Use this command to actually run the merge_tops_burst workflow
get_ipython().system('python -m hyp3_isce2 ++process merge_tops_bursts $data_dir')


# ## 3. Display the merged product
# 
# We've successfully run the command, now let's look at the results!

# In[ ]:


get_ipython().run_line_magic('matplotlib', 'inline')
import numpy as np
from osgeo import gdal
import matplotlib
import matplotlib.pyplot as plt

tifs = [f for f in work_dir.glob("*/*.tif")]
unw_file = [f for f in tifs if 'unw_phase' in f.name][0]
wrapped_file = [f for f in tifs if 'wrapped_phase' in f.name][0]
corr_file = [f for f in tifs if 'corr' in f.name][0]
desired_tifs = [wrapped_file, unw_file, corr_file]

f, axs = plt.subplots(len(desired_tifs), figsize=(6,10))
for i, tif in enumerate(desired_tifs): 
    ds = gdal.Open(str(tif))
    merged_bursts = np.ma.masked_equal(ds.GetRasterBand(1).ReadAsArray(), 0)
    ds = None
        
    axs[i].imshow(merged_bursts)
    axs[i].set_title(tif.name)

    plt.setp(plt.gcf().get_axes(), xticks=[], yticks=[])
    plt.tight_layout()


# In[ ]: