topsStack_parallel.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<font size=\"1\">Copyright 2021, by the California Institute of Technology. ALL RIGHTS RESERVED. United States Government sponsorship acknowledged. Any commercial use must be negotiated with the Office of Technology Transfer at the California Institute of Technology.</font>\n",
    "    \n",
    "<font size=\"1\">This software may be subject to U.S. export control laws and regulations. By accepting this document, the user agrees to comply with all applicable U.S. export laws and regulations. User has the responsibility to obtain export licenses, or other export authority as may be required, before exporting such information to foreign countries or providing access to foreign persons.<font>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Sentinel-1 TOPS stack processor\n",
    "The detailed algorithm for stack processing of TOPS data can be find here:\n",
    "\n",
    "+ Fattahi, H., P. Agram, and M. Simons (2016), A Network-Based Enhanced Spectral Diversity Approach for TOPS Time-Series Analysis, IEEE Transactions on Geoscience and Remote Sensing, 55(2), 777-786, doi:[10.1109/TGRS.2016.2614925](https://ieeexplore.ieee.org/abstract/document/7637021).\n",
    "\n",
    "-----------------------------------\n",
    "\n",
    "The scripts provides support for Sentinel-1 TOPS stack processing. Currently supported workflows include a coregistered stack of SLC, interferograms, offsets, and coherence. \n",
    "\n",
    "#### Be sure [default]  the ~/.netrc file has been changed to include valid credentials for urs.earthdata.nasa.gov."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_current_time():\n",
    "    import datetime\n",
    "    return datetime.datetime.now()\n",
    "run_start_time = get_current_time()\n",
    "print(\"PGE run start time : {}\".format(run_start_time))\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import os\n",
    "import json\n",
    "from math import floor, ceil\n",
    "import json\n",
    "import re\n",
    "import osaka\n",
    "import osaka.main\n",
    "from builtins import str\n",
    "import os, sys, re, json, logging, traceback, requests, argparse\n",
    "from datetime import datetime\n",
    "from pprint import pformat\n",
    "from requests.packages.urllib3.exceptions import (InsecureRequestWarning,\n",
    "                                                  InsecurePlatformWarning)\n",
    "try: from html.parser import HTMLParser\n",
    "except: from html.parser import HTMLParser\n",
    "        \n",
    "requests.packages.urllib3.disable_warnings(InsecureRequestWarning)\n",
    "requests.packages.urllib3.disable_warnings(InsecurePlatformWarning)\n",
    "PROCESSING_START=datetime.utcnow().isoformat()\n",
    "print(PROCESSING_START)\n",
    "\n",
    "# this block makes sure the directory set-up/change is only done once and relative to the notebook's directory\\n\",\n",
    "try:\n",
    "    start_dir\n",
    "except NameError:\n",
    "    start_dir = os.getcwd()\n",
    "    output_dir = os.path.join(start_dir, 'notebook_output/topsStack_parallel')\n",
    "    os.makedirs(output_dir, exist_ok=True)\n",
    "\n",
    "os.chdir(output_dir)\n",
    "    \n",
    "PGE_BASE=os.getcwd()\n",
    "SLC_PATH = \"zip\"\n",
    "print(PGE_BASE)\n",
    "\n",
    "SLC_RE = re.compile(r'(?P<mission>S1\\w)_IW_SLC__.*?' +\n",
    "                    r'_(?P<start_year>\\d{4})(?P<start_month>\\d{2})(?P<start_day>\\d{2})' +\n",
    "                    r'T(?P<start_hour>\\d{2})(?P<start_min>\\d{2})(?P<start_sec>\\d{2})' +\n",
    "                    r'_(?P<end_year>\\d{4})(?P<end_month>\\d{2})(?P<end_day>\\d{2})' +\n",
    "                    r'T(?P<end_hour>\\d{2})(?P<end_min>\\d{2})(?P<end_sec>\\d{2})_.*$')\n",
    "ISCE_HOME=\"/opt/isce2/isce\"\n",
    "\n",
    "QC_SERVER = 'https://qc.sentinel1.eo.esa.int/'\n",
    "DATA_SERVER = 'http://aux.sentinel1.eo.esa.int/'\n",
    "\n",
    "ORBITMAP = [('precise','aux_poeorb', 100),\n",
    "            ('restituted','aux_resorb', 100)]\n",
    "\n",
    "OPER_RE = re.compile(r'S1\\w_OPER_AUX_(?P<type>\\w+)_OPOD_(?P<yr>\\d{4})(?P<mo>\\d{2})(?P<dy>\\d{2})')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "tags": [
     "parameters"
    ]
   },
   "outputs": [],
   "source": [
    "from typing import List\n",
    "ISCE_HOME=\"/opt/isce2/isce\"\n",
    "min_lat = 34.6002832\n",
    "max_lat = 34.6502392\n",
    "min_lon = -79.0801608\n",
    "max_lon = -78.9705888\n",
    "master_date = \"\"\n",
    "localize_slcs: List[str] = [\"S1A_IW_SLC__1SDV_20150315T231319_20150315T231349_005049_006569_0664\",\n",
    "                      \"S1A_IW_SLC__1SDV_20150818T231326_20150818T231356_007324_00A0E0_93D5\",\n",
    "                      \"S1A_IW_SLC__1SDV_20150830T231332_20150830T231402_007499_00A5A9_02B3\",\n",
    "                      \"S1A_IW_SLC__1SDV_20160414T231321_20160414T231350_010824_01030B_F02B\",\n",
    "                      \"S1A_IW_SLC__1SDV_20160414T231348_20160414T231416_010824_01030B_9FA9\"\n",
    "                     ] \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ctx = {}\n",
    "ctx[\"min_lat\"] = min_lat\n",
    "ctx[\"max_lat\"] = max_lat\n",
    "ctx[\"min_lon\"] = min_lon\n",
    "ctx[\"max_lon\"] = max_lon\n",
    "ctx[\"master_date\"]=\"\"\n",
    "ctx[\"localize_slcs\"]= localize_slcs\n",
    "\n",
    "print(json.dumps(ctx, indent=4))\n",
    "wd = os.getcwd()\n",
    "global runtime_dict\n",
    "runtime_dict = {}\n",
    "\n",
    "with open('_stdout.txt', 'w') as f:\n",
    "    f.write(\"Output File\")\n",
    "with open('_context.json', 'w') as outfile:\n",
    "    json.dump(ctx, outfile)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import requests\n",
    "import re\n",
    "import os\n",
    "from datetime import datetime, timedelta\n",
    "import urllib.parse\n",
    "\n",
    "def session_get(session, url):\n",
    "    return session.get(url, verify=False)\n",
    "\n",
    "datefmt = \"%Y%m%dT%H%M%S\"\n",
    "queryfmt = \"%Y-%m-%d\"\n",
    "server = 'https://scihub.copernicus.eu/gnss/'\n",
    "auth_netloc_fmt = \"{}:{}@{}\"\n",
    "url_tpt = 'search?q=( beginPosition:[{0}T00:00:00.000Z TO {1}T23:59:59.999Z] AND endPosition:[{0}T00:00:00.000Z TO {1}T23:59:59.999Z] ) AND ( (platformname:Sentinel-1 AND filename:{2}_* AND producttype:{3}))&start=0&rows=100'\n",
    "credentials = ('gnssguest','gnssguest')\n",
    "\n",
    "class MyHTMLParser(HTMLParser):\n",
    "\n",
    "    def __init__(self,url):\n",
    "        HTMLParser.__init__(self)\n",
    "        self.fileList = []\n",
    "        self._url = url\n",
    "\n",
    "    def handle_starttag(self, tag, attrs):\n",
    "        for name, val in attrs:\n",
    "            if name == 'href':\n",
    "                if val.startswith(\"https://scihub.copernicus.eu/gnss/odata\") and val.endswith(\")/\"):\n",
    "                    pass\n",
    "                else:\n",
    "                    downloadLink = val.strip()\n",
    "                    downloadLink = downloadLink.split(\"/Products('Quicklook')\")\n",
    "                    downloadLink = downloadLink[0] + downloadLink[-1]\n",
    "                    self._url = downloadLink\n",
    "\n",
    "    def handle_data(self, data):\n",
    "        if data.startswith(\"S1\") and data.endswith(\".EOF\"):\n",
    "            self.fileList.append((self._url, data.strip()))\n",
    "\n",
    "\n",
    "def fileToRange(fname):\n",
    "    '''\n",
    "    Derive datetime range from orbit file name.\n",
    "    '''\n",
    "\n",
    "    fields = os.path.basename(fname).split('_')\n",
    "    start = datetime.strptime(fields[-2][1:16], datefmt)\n",
    "    stop = datetime.strptime(fields[-1][:15], datefmt)\n",
    "    mission = fields[0]\n",
    "\n",
    "    return (start, stop, mission)\n",
    "\n",
    "\n",
    "def get_download_orbit_dict(slc_start_dt, slc_end_dt, sat_name):\n",
    "\n",
    "    found = False\n",
    "    delta = timedelta(days=1)\n",
    "    timebef = (slc_end_dt - delta).strftime(queryfmt)\n",
    "    timeaft = (slc_end_dt + delta).strftime(queryfmt)\n",
    "    match = None\n",
    "    matchFileName = None\n",
    "    \n",
    "    session = requests.Session()\n",
    "\n",
    "    for fidelity in ('AUX_POEORB', 'AUX_RESORB'):\n",
    "        url = server + url_tpt.format(timebef, timeaft, sat_name, fidelity)\n",
    "        # print(f\"url: {url}\")\n",
    "\n",
    "        try:\n",
    "            r = session.get(url, verify=True, auth=credentials)\n",
    "            r.raise_for_status()\n",
    "            parser = MyHTMLParser(url)\n",
    "            parser.feed(r.text)\n",
    "\n",
    "            for resulturl, result in parser.fileList:\n",
    "                # print('Results: {} : {}'.format(resulturl, result))\n",
    "                tbef, taft, mission = fileToRange(os.path.basename(result))\n",
    "\n",
    "                if (tbef <= slc_start_dt) and (taft >= slc_end_dt):\n",
    "                    matchFileName = result\n",
    "                    parse_url = urllib.parse.urlsplit(resulturl)\n",
    "                    # Add credentials for osaka\n",
    "                    new_netloc = auth_netloc_fmt.format(credentials[0],credentials[1],parse_url.netloc)\n",
    "                    match = urllib.parse.urlunsplit(parse_url._replace(netloc=new_netloc))\n",
    "                else:\n",
    "                    print(\"no match.\")\n",
    "            \n",
    "            if match is not None:\n",
    "                found = True\n",
    "        except Exception as e:\n",
    "            print(\"Exception {}\".format(e))\n",
    "        \n",
    "        if found:\n",
    "            break\n",
    "\n",
    "    # print(\"returning {} : {}\".format(matchFileName, match))\n",
    "    return matchFileName, match\n",
    "\n",
    "\n",
    "def get_orbit_files():\n",
    "    import json\n",
    "    import os\n",
    "    import osaka\n",
    "    import osaka.main\n",
    "    orbit_dates = []\n",
    "        \n",
    "    for slc in ctx[\"localize_slcs\"]:\n",
    "        match = SLC_RE.search(slc)\n",
    "        if not match:\n",
    "            raise RuntimeError(\"Failed to recognize SLC ID %s.\" %slc)\n",
    "        mission = match.group('mission')\n",
    "        slc_start_dt_str = \"{}-{}-{}T{}:{}:{}\".format(match.group('start_year'), \n",
    "                                                      match.group('start_month'),\n",
    "                                                      match.group('start_day'),\n",
    "                                                      match.group('start_hour'),\n",
    "                                                      match.group('start_min'),\n",
    "                                                      match.group('start_sec'))\n",
    "        slc_start_dt = datetime.strptime(slc_start_dt_str, \"%Y-%m-%dT%H:%M:%S\")\n",
    "        slc_end_dt_str = \"{}-{}-{}T{}:{}:{}\".format(match.group('end_year'), \n",
    "                                                      match.group('end_month'),\n",
    "                                                      match.group('end_day'),\n",
    "                                                      match.group('end_hour'),\n",
    "                                                      match.group('end_min'),\n",
    "                                                      match.group('end_sec'))\n",
    "        slc_end_dt = datetime.strptime(slc_end_dt_str, \"%Y-%m-%dT%H:%M:%S\")\n",
    "        day_dt_str = slc_start_dt.strftime('%Y-%m-%d')\n",
    "        if day_dt_str not in orbit_dates:\n",
    "            orbit_dates.append(day_dt_str)\n",
    "            orbit_file_name, orbit_file_uri = get_download_orbit_dict(slc_start_dt, slc_end_dt, mission)\n",
    "            if orbit_file_uri is not None:\n",
    "                directory = os.path.join(wd, \"orbits\")\n",
    "                if not os.path.exists(directory):\n",
    "                    os.makedirs(directory)\n",
    "                    print(f\"created {directory}\")\n",
    "                osaka.main.get(orbit_file_uri, os.path.join(directory, orbit_file_name))\n",
    "                print(\"Downloaded orbit file {}\".format(orbit_file_name))\n",
    "            else:\n",
    "                print('No orbit files found for slc: {}'.format(slc))\n",
    "\n",
    "\n",
    "def download_slc(slc_id, path):\n",
    "    url = \"https://datapool.asf.alaska.edu/SLC/SA/{}.zip\".format(slc_id)\n",
    "    print(\"Downloading {} : {} to {}\".format(slc_id, url, path))\n",
    "    \n",
    "    if not os.path.exists(path):\n",
    "        os.makedirs(path)\n",
    "    osaka.main.get(url, path)\n",
    "\n",
    "def run_cmd_output(cmd):\n",
    "    from subprocess import check_output, CalledProcessError\n",
    "    cmd_line = \" \".join(cmd)\n",
    "    print(\"Calling: {}\".format(cmd_line))\n",
    "    output = check_output(cmd_line, shell=True)\n",
    "    return output\n",
    "\n",
    "def run_cmd(cmd):\n",
    "    !source /opt/isce2/isce_env.sh\n",
    "    !export PATH=/opt/conda/bin/:/opt/isce2/src/isce2/contrib/stack/topsStack/:${PATH}\n",
    "    !export PYTHONPATH=/opt/isce2/src/isce2/contrib/stack/topsStack/:${PYTHONPATH}\n",
    "\n",
    "    import subprocess\n",
    "    from subprocess import check_call, CalledProcessError\n",
    "    import sys\n",
    "    print(cmd)\n",
    "    cmd_line = \" \".join(cmd)\n",
    "    print(\"Calling : {}\".format(cmd_line))\n",
    "    p = subprocess.Popen(cmd_line, shell=True,stdout=subprocess.PIPE,stderr=subprocess.STDOUT)\n",
    "    while True: \n",
    "        line = p.stdout.readline()\n",
    "        if not line:\n",
    "            break\n",
    "        print(line.strip())\n",
    "        sys.stdout.flush()\n",
    "        \n",
    "def get_minimum_bounding_rectangle():\n",
    "    from math import floor, ceil\n",
    "    cwd = os.getcwd()\n",
    "    slc_ids = [x for x in os.listdir('.') if os.path.isdir(x) and '_SLC__' in x]\n",
    "\n",
    "    all_lats = []\n",
    "    all_lons = []\n",
    "    for slc in slc_ids:\n",
    "        slc_met_json = slc + '.met.json'\n",
    "\n",
    "        with open(os.path.join(cwd, slc, slc_met_json), 'r') as f:\n",
    "            data = json.load(f)\n",
    "            bbox = data['bbox']\n",
    "            for coord in bbox:\n",
    "                all_lats.append(coord[0])\n",
    "                all_lons.append(coord[1])\n",
    "\n",
    "    min_lat = min(all_lats) + 0.2\n",
    "    max_lat = max(all_lats) - 0.1\n",
    "    min_lon = min(all_lons)\n",
    "    max_lon = max(all_lons)\n",
    "\n",
    "    min_lat_lo = floor(min_lat)\n",
    "    max_lat_hi = ceil(max_lat)\n",
    "    min_lon_lo = floor(min_lon)\n",
    "    max_lon_hi = ceil(max_lon)\n",
    "\n",
    "    return min_lat, max_lat, min_lon, max_lon, min_lat_lo, max_lat_hi, min_lon_lo, max_lon_hi\n",
    "\n",
    "def get_user_input_bbox(ctx_file):\n",
    "    \"\"\"\n",
    "    :param ctx_file: dictionary from cxt file\n",
    "    :return: void\n",
    "    \"\"\"\n",
    "    from math import floor, ceil\n",
    "    min_lat = ctx_file['min_lat']\n",
    "    max_lat = ctx_file['max_lat']\n",
    "    min_lon = ctx_file['min_lon']\n",
    "    max_lon = ctx_file['max_lon']\n",
    "\n",
    "    min_lat_lo = floor(min_lat)\n",
    "    max_lat_hi = ceil(max_lat)\n",
    "    min_lon_lo = floor(min_lon)\n",
    "    max_lon_hi = ceil(max_lon)\n",
    "\n",
    "    return min_lat, max_lat, min_lon, max_lon, min_lat_lo, max_lat_hi, min_lon_lo, max_lon_hi\n",
    "\n",
    "def get_master_date(ctx):\n",
    "    master_date = ctx.get('master_date', \"\")\n",
    "    return master_date\n",
    "\n",
    "def get_bbox(ctx):\n",
    "    # min_lat, max_lat, min_lon, max_lon = ctx['region_of_interest']\n",
    "\n",
    "    if ctx['min_lat'] != \"\" or ctx['max_lat'] != \"\" or ctx['min_lon'] != \"\" or ctx['max_lon'] != \"\":\n",
    "        # if any values are present in _context.json we can assume user put them in manually\n",
    "        bbox_data = get_user_input_bbox(ctx)\n",
    "    else:\n",
    "        # if user did not define ANY lat lons\n",
    "        bbox_data = get_minimum_bounding_rectangle()\n",
    "\n",
    "    return bbox_data\n",
    "\n",
    "def download_dem():\n",
    "    dem_cmd = [\n",
    "        \"{}/applications/dem.py\".format(ISCE_HOME), \"-a\",\n",
    "        \"stitch\", \"-b\", \"{} {} {} {}\".format(MINLAT_LO, MAXLAT_HI, MINLON_LO, MAXLON_HI),\n",
    "        \"-r\", \"-s\", \"1\", \"-f\", \"-c\", \"|\", \"tee\", \"dem.txt\"\n",
    "        #\"-n\", dem_user, \"-w\", dem_pass,\"-u\", dem_url\n",
    "    ]\n",
    "    run_cmd(dem_cmd)\n",
    "    \n",
    "def get_config_file(run_file, kward):\n",
    "    config_file = None\n",
    "    with open(run_file, 'r') as fp:\n",
    "        line = fp.readline()\n",
    "        while line:\n",
    "            if kward in line:\n",
    "                config_file = line.split(' ')[-1]\n",
    "                break\n",
    "            line = fp.readline()\n",
    "    return config_file.strip()\n",
    "\n",
    "def get_config_val(config_file, p_name):\n",
    "    \n",
    "    print(\"config_file : {}, p_name : {}\".format(config_file, p_name))\n",
    "    val = None\n",
    "    with open(config_file, 'r') as fp:\n",
    "        line = fp.readline()\n",
    "        \n",
    "        while line:\n",
    "            if p_name in line:\n",
    "                val = line.split(':')[-1].strip()\n",
    "                break\n",
    "            line = fp.readline()\n",
    "    return val\n",
    "\n",
    "def get_date_from_str(test_str):\n",
    "    import re\n",
    "    temp = re.findall(r'\\d+', test_str) \n",
    "    res = list(map(int, temp))\n",
    "    return str(res[0])\n",
    "\n",
    "\n",
    "def stackSlcDn_run2_5():\n",
    "    run1File = \"./run_files/run_01_unpack_topo_reference\"\n",
    "    run2File = \"./run_files/run_02_unpack_secondary_slc\"\n",
    "    runFile = \"./run_files/run_02.5_slc_noise_calibration\"\n",
    "\n",
    "    run_cmd(\"rm -f {}\".format(runFile).split(' '))\n",
    "    run_cmd(\"touch {}\".format(runFile).split(' '))\n",
    "    reference_dir = os.path.join(wd, \"reference\")\n",
    "    secondary_dir = os.path.join(wd, \"secondary\")\n",
    "    refConfig = get_config_file(run1File, \"reference\")\n",
    "    print(refConfig)\n",
    "    \n",
    "    ref_zip = get_config_val(refConfig, 'dirname')\n",
    "    print(ref_zip)\n",
    "    ref_swath = get_config_val(refConfig, 'swaths')\n",
    "    print(ref_swath)\n",
    "    \n",
    "    with open(runFile, 'w') as fw:\n",
    "        cmd = \"read_calibration_slc.py -zip {} -ext {} -od {} -o -t noise -n '{}'\\n\".format(ref_zip, bbox, reference_dir, ref_swath)\n",
    "        fw.write(cmd)\n",
    "        with open(run2File, 'r') as fp:\n",
    "            line = fp.readline()\n",
    "            cmds = []\n",
    "            while line:\n",
    "                if 'secondary' in line:\n",
    "                    print(line)\n",
    "                    secConfig = line.split(' ')[-1].strip()\n",
    "                    sec_date = get_date_from_str(os.path.basename(secConfig))\n",
    "                    print(\"{} : {}\".format(sec_date, os.path.basename(secConfig)))\n",
    "                    sec_zip = get_config_val(secConfig, 'dirname')\n",
    "                    sec_swath = get_config_val(secConfig, 'swaths')\n",
    "                    odir=os.path.join(secondary_dir, sec_date )\n",
    "                    print(\"{} : {} :{}\".format(sec_zip, sec_swath, odir ))\n",
    "                    cmd = \"read_calibration_slc.py -zip {} -ext {} -od {} -o -t noise -n '{}'\\n\".format(sec_zip, bbox, odir, sec_swath)\n",
    "                    fw.write(cmd)\n",
    "            \n",
    "                line = fp.readline()\n",
    "\n",
    "\n",
    "    \n",
    "import os\n",
    "import sys\n",
    "import shutil\n",
    "import json\n",
    "import re\n",
    "\n",
    "from datetime import datetime\n",
    "from osgeo import ogr, osr\n",
    "\n",
    "def get_union_polygon_from_bbox(env):\n",
    "    print(env)\n",
    "    if not isinstance(env, list):\n",
    "        env = env.split()\n",
    "    env = [float(i) for i in env] \n",
    "    coords = [\n",
    "        [ env[3], env[0] ],\n",
    "        [ env[3], env[1] ],\n",
    "        [ env[2], env[1] ],\n",
    "        [ env[2], env[0] ],\n",
    "        [ env[3], env[0] ]\n",
    "    ]\n",
    "    return {\n",
    "        \"type\": \"polygon\",\n",
    "        \"coordinates\": [ coords ]\n",
    "    }\n",
    "\n",
    "# copied from stitch_ifgs.get_union_polygon()\n",
    "def get_union_polygon(ds_files):\n",
    "    \"\"\"\n",
    "    Get GeoJSON polygon of union of IFGs.\n",
    "    :param ds_files: list of .dataset.json files, which have the 'location' key\n",
    "    :return: geojson of merged bbox\n",
    "    \"\"\"\n",
    "\n",
    "    geom_union = None\n",
    "    for ds_file in ds_files:\n",
    "        f = open(ds_file)\n",
    "        ds = json.load(f)\n",
    "        geom = ogr.CreateGeometryFromJson(json.dumps(ds['location'], indent=2, sort_keys=True))\n",
    "        if geom_union is None:\n",
    "            geom_union = geom\n",
    "        else:\n",
    "            geom_union = geom_union.Union(geom)\n",
    "    return json.loads(geom_union.ExportToJson()), geom_union.GetEnvelope()\n",
    "\n",
    "\n",
    "def get_dataset_met_json_files(cxt):\n",
    "    \"\"\"\n",
    "    returns 2 lists: file paths for dataset.json files and met.json files\n",
    "    :param cxt: json from _context.json\n",
    "    :return: list[str], list[str]\n",
    "    \"\"\"\n",
    "    pwd = os.getcwd()\n",
    "    localize_urls = cxt['localize_slcs']\n",
    "\n",
    "    met_files, ds_files = [], []\n",
    "    for localize_url in localize_urls:\n",
    "        local_path = localize_url['local_path']\n",
    "        slc_id = local_path.split('/')[0]\n",
    "        slc_path = os.path.join(pwd, slc_id, slc_id)\n",
    "\n",
    "        ds_files.append(slc_path + '.dataset.json')\n",
    "        met_files.append(slc_path + '.met.json')\n",
    "    return ds_files, met_files\n",
    "\n",
    "\n",
    "def get_scenes(cxt):\n",
    "    \"\"\"\n",
    "    gets all SLC scenes for the stack\n",
    "    :param cxt: contents for _context.json\n",
    "    :return: list of scenes\n",
    "    \"\"\"\n",
    "    localize_urls = cxt['localize_slcs']\n",
    "    all_scenes = set()\n",
    "    for localize_url in localize_urls:\n",
    "        local_path = localize_url['local_path']\n",
    "        slc_id = local_path.split('/')[0]\n",
    "        all_scenes.add(slc_id)\n",
    "    return sorted(list(all_scenes))\n",
    "\n",
    "\n",
    "def get_min_max_timestamps(scenes_ls):\n",
    "    \"\"\"\n",
    "    returns the min timestamp and max timestamp of the stack\n",
    "    :param scenes_ls: list[str] all slc scenes in stack\n",
    "    :return: (str, str) 2 timestamp strings, ex. 20190518T161611\n",
    "    \"\"\"\n",
    "    timestamps = set()\n",
    "\n",
    "    regex_pattern = r'(\\d{8}T\\d{6}).(\\d{8}T\\d{6})'\n",
    "    for scene in scenes_ls:\n",
    "        matches = re.search(regex_pattern, scene)\n",
    "        if not matches:\n",
    "            raise Exception(\"regex %s was unable to match with SLC id %s\" % (regex_pattern, scene))\n",
    "\n",
    "        slc_timestamps = (matches.group(1), matches.group(2))\n",
    "        timestamps = timestamps.union(slc_timestamps)\n",
    "\n",
    "    #min_timestamp = \"{}Z\".format(datetime.strptime(min(timestamps), \"%Y%m%dT%H%M%S\").isoformat())\n",
    "    #max_timestamp = \"{}Z\".format(datetime.strptime(max(timestamps), \"%Y%m%dT%H%M%S\").isoformat())\n",
    "    \n",
    "    return min(timestamps), max(timestamps)\n",
    "\n",
    "def get_min_max_times(scenes_ls):\n",
    "    \"\"\"\n",
    "    returns the min timestamp and max timestamp of the stack\n",
    "    :param scenes_ls: list[str] all slc scenes in stack\n",
    "    :return: (str, str) 2 timestamp strings, ex. 20190518T161611\n",
    "    \"\"\"\n",
    "    timestamps = set()\n",
    "\n",
    "    regex_pattern = r'(\\d{8}T\\d{6}).(\\d{8}T\\d{6})'\n",
    "    for scene in scenes_ls:\n",
    "        matches = re.search(regex_pattern, scene)\n",
    "        if not matches:\n",
    "            raise Exception(\"regex %s was unable to match with SLC id %s\" % (regex_pattern, scene))\n",
    "\n",
    "        slc_timestamps = (matches.group(1), matches.group(2))\n",
    "        timestamps = timestamps.union(slc_timestamps)\n",
    "\n",
    "    #min_timestamp = min(timestamps)\n",
    "    #max_timestamp = max(timestamps)\n",
    "    min_timestamp = \"{}\".format(datetime.strptime(min(timestamps), \"%Y%m%dT%H%M%S\").isoformat())\n",
    "    max_timestamp = \"{}\".format(datetime.strptime(max(timestamps), \"%Y%m%dT%H%M%S\").isoformat())\n",
    "    \n",
    "    return min_timestamp, max_timestamp\n",
    "\n",
    "\n",
    "def create_list_from_keys_json_file(json_files, *args):\n",
    "    \"\"\"\n",
    "    gets all key values in each .json file and returns a sorted array of values\n",
    "    :param json_files: list[str]\n",
    "    :return: list[]\n",
    "    \"\"\"\n",
    "    values = set()\n",
    "    for json_file in json_files:\n",
    "        if os.path.isfile(json_file): \n",
    "            f = open(json_file)\n",
    "            data = json.load(f)\n",
    "            for arg in args:\n",
    "                value = data[arg]\n",
    "                values.add(value)\n",
    "    return sorted(list(values))\n",
    "\n",
    "\n",
    "def camelcase_to_underscore(name):\n",
    "    s1 = re.sub('(.)([A-Z][a-z]+)', r'\\1_\\2', name)\n",
    "    return re.sub('([a-z0-9])([A-Z])', r'\\1_\\2', s1).lower()\n",
    "\n",
    "\n",
    "def get_key_and_convert_to_underscore(json_file_paths, key):\n",
    "    \"\"\"\n",
    "    read through all the json files in file paths, get the first occurrence of key and convert it to underscore\n",
    "    :param json_file_paths: list[str]\n",
    "    :param key: str\n",
    "    :return: key and value\n",
    "    \"\"\"\n",
    "    for json_file in json_file_paths:\n",
    "        if os.path.isfile(json_file): \n",
    "            f = open(json_file)\n",
    "            data = json.load(f)\n",
    "            if key in data.keys():\n",
    "                underscore_key = camelcase_to_underscore(key)\n",
    "                return underscore_key, data[key]\n",
    "    return None, None\n",
    "\n",
    "\n",
    "def generate_dataset_json_data(ctx, version):\n",
    "    \"\"\"\n",
    "    :param cxt: _context.json file\n",
    "    :param dataset_json_files: list[str] all file paths of SLC's .dataset.json files\n",
    "    :param version: str: version, ex. v1.0\n",
    "    :return: dict\n",
    "    \"\"\"\n",
    "    dataset_json_data = dict()\n",
    "    dataset_json_data['version'] = version\n",
    "\n",
    "    \n",
    "    try:      \n",
    "        dataset_json_data['starttime'], dataset_json_data['endtime']  = get_min_max_times(ctx[\"localize_slcs\"])      \n",
    "    except Exception as err:\n",
    "        print(str(err))\n",
    "  \n",
    "    try:\n",
    "        dataset_json_data['location'] = get_union_polygon_from_bbox(bbox.split(' '))\n",
    "    except Exception as err:\n",
    "        print(str(err))\n",
    "\n",
    "    return dataset_json_data\n",
    "\n",
    "\n",
    "def generate_met_json_data(ctx, bbox, version):\n",
    "    \"\"\"\n",
    "    :param cxt: _context.json file\n",
    "    :param met_json_file_paths: list[str] all file paths of SLC's .met.json files\n",
    "    :param dataset_json_files: list[str] all file paths of SLC's .dataset.json files\n",
    "    :param version: str: version, ex. v1.0\n",
    "    :return: dict\n",
    "    \"\"\"\n",
    "    met_json_data = {\n",
    "        'processing_start': '{}'.format(PROCESSING_START),\n",
    "        'processing_stop': '{}'.format(datetime.utcnow().isoformat()),\n",
    "        'version': version\n",
    "    }\n",
    "\n",
    "    # generating bbox\n",
    "    geojson = get_union_polygon_from_bbox(bbox)\n",
    "    coordinates = geojson['coordinates'][0]\n",
    "    for coordinate in coordinates:\n",
    "        coordinate[0], coordinate[1] = coordinate[1], coordinate[0]\n",
    "    \n",
    "    met_json_data['bbox'] = bbox\n",
    "    \n",
    "    # list of SLC scenes\n",
    "   \n",
    "    met_json_data['scenes'] = ctx[\"localize_slcs\"]\n",
    "    met_json_data['scene_count'] = len(ctx[\"localize_slcs\"])\n",
    "\n",
    "    # getting timestamps\n",
    "\n",
    "    met_json_data['sensing_start'], met_json_data['sensing_stop']  = get_min_max_times(ctx[\"localize_slcs\"])\n",
    "\n",
    "    # additional information\n",
    "    met_json_data['dataset_type'] = 'topsStack_hamsar'\n",
    "\n",
    "    return met_json_data\n",
    "\n",
    "\n",
    "def read_context():\n",
    "    with open('_context.json', 'r') as f:\n",
    "        cxt = json.load(f)\n",
    "        return cxt\n",
    "\n",
    "\n",
    "def create_dataset(bbox):\n",
    "    VERSION = 'v1.0'\n",
    "    DATASET_NAMING_TEMPLATE = 'coregistered_slcs-{min_timestamp}-{max_timestamp}'\n",
    "    PWD = os.getcwd()\n",
    "\n",
    "    # creating list of all SLC .dataset.json and .met.json files\n",
    "    ctx = read_context()\n",
    "    #dataset_json_files, met_json_files = get_dataset_met_json_files(context_json)\n",
    "\n",
    "    # getting list of SLC scenes and extracting min max timestamp\n",
    "    slc_scenes = ctx['localize_slcs']\n",
    "    min_timestamp, max_timestamp = get_min_max_timestamps(slc_scenes)\n",
    "\n",
    "    # creatin dataset directory\n",
    "    dataset_name = DATASET_NAMING_TEMPLATE.format(min_timestamp=min_timestamp, max_timestamp=max_timestamp)\n",
    "    if not os.path.exists(dataset_name):\n",
    "        os.mkdir(dataset_name)\n",
    "\n",
    "    # move merged/ master/ slaves/ directory to dataset directory\n",
    "    move_directories = ['merged', 'reference', 'secondarys']\n",
    "    for directory in move_directories:\n",
    "        shutil.move(directory, dataset_name)\n",
    "\n",
    "    # move _stdout.txt log file to dataset\n",
    "    shutil.copyfile('_stdout.txt', os.path.join(dataset_name, '_stdout.txt'))\n",
    "\n",
    "    # generate .dataset.json data\n",
    "    dataset_json_data = {}\n",
    "    try:\n",
    "        dataset_json_data = generate_dataset_json_data(ctx, VERSION)\n",
    "    except Exception as err:\n",
    "        print(str(err))\n",
    "    dataset_json_data['label'] = dataset_name\n",
    "    print(json.dumps(dataset_json_data, indent=2))\n",
    "\n",
    "    # generate .met.json data\n",
    "    met_json_data = generate_met_json_data(ctx, bbox, VERSION)\n",
    "    print(json.dumps(met_json_data, indent=2))\n",
    "\n",
    "    # writing .dataset.json to file\n",
    "    dataset_json_filename = os.path.join(PWD, dataset_name, dataset_name + '.dataset.json')\n",
    "    with open(dataset_json_filename, 'w') as f:\n",
    "        json.dump(dataset_json_data, f, indent=2)\n",
    "\n",
    "    # writing .met.json to file\n",
    "    met_json_filename = os.path.join(PWD, dataset_name, dataset_name + '.met.json')\n",
    "    with open(met_json_filename, 'w') as f:\n",
    "        json.dump(met_json_data, f, indent=2)\n",
    "        \n",
    "        \n",
    "        \n",
    "def download_slcs(path):\n",
    "    \n",
    "    for slc in ctx[\"localize_slcs\"]:\n",
    "        download_slc(slc, path)       "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import xml.etree.ElementTree as ET\n",
    "import argparse\n",
    "from argparse import RawTextHelpFormatter\n",
    "import zipfile\n",
    "import fnmatch\n",
    "import re\n",
    "import os\n",
    "import sys        \n",
    "import isce\n",
    "import iscesys\n",
    "import isceobj\n",
    "import isceobj.Sensor.TOPS as TOPS\n",
    "import isceobj.Sensor.TOPS.BurstSLC as BurstSLC\n",
    "from isceobj.Image import createImage\n",
    "import glob\n",
    "from osgeo import gdal\n",
    "import numpy as np\n",
    "from scipy.interpolate import LinearNDInterpolator as interpnd\n",
    "from iscesys.Parsers import XmlParser\n",
    "\n",
    "\n",
    "\n",
    "def locateCaliFile(slc,type,polid='vv'):\n",
    "    swathid = 's1?-iw%d'%(slc.swathNumber)\n",
    "#    polid = slc.polarization  #default is already set to 'vv'\n",
    "    print('Using data polarization ', polid)\n",
    "    match = None\n",
    "    for dirname in slc.safe:\n",
    "        match = None\n",
    "        \n",
    "        if dirname.endswith('.zip'):\n",
    "            zf = zipfile.ZipFile(dirname, 'r')\n",
    "            if type == 'radio':\n",
    "                pattern = os.path.join('*SAFE','annotation','calibration','calibration-') + swathid + '-slc-' + polid + '*.xml'\n",
    "                match = fnmatch.filter(zf.namelist(), pattern)\n",
    "                if (len(match) == 0):\n",
    "                    raise Exception('No radiometric calibration file found in zip file: {0}'.format(dirname))\n",
    "                slc.radioCali.append('/vsizip/'+os.path.join(dirname, match[0]) )\n",
    "                print('Found radiometric calibration files: ', slc.radioCali)\n",
    "            elif type == 'noise':\n",
    "                pattern = os.path.join('*SAFE','annotation','calibration','noise-') + swathid + '-slc-' + polid + '*.xml'\n",
    "                zf = zipfile.ZipFile(dirname, 'r')\n",
    "                match = fnmatch.filter(zf.namelist(), pattern)\n",
    "                if (len(match) == 0):\n",
    "                    raise Exception('No noise calibration file found in zip file: {0}'.format(dirname))\n",
    "                slc.noiseCali.append('/vsizip/'+os.path.join(dirname, match[0]) )\n",
    "                print('Found noise calibration files: ', slc.noiseCali)\n",
    "            zf.close()\n",
    "        \n",
    "        else:\n",
    "            if type == 'radio':\n",
    "                pattern = os.path.join('annotation','calibration','calibration-') + swathid + '-slc-' + polid + '*.xml'\n",
    "                match = glob.glob( os.path.join(dirname, pattern))\n",
    "                if (len(match) == 0):\n",
    "                    raise Exception('No radiometric calibration file found in {0}'.format(dirname))\n",
    "                slc.radioCali.append(match[0])\n",
    "                print('Found radiometric calibration files: ', slc.radioCali)\n",
    "            elif type == 'noise':\n",
    "                pattern = os.path.join('annotation','calibration','noise-') + swathid + '-slc-' + polid + '*.xml'\n",
    "                match = glob.glob( os.path.join(dirname, pattern))\n",
    "                if (len(match) == 0):\n",
    "                    raise Exception('No noise calibration file found in {0}'.format(dirname))\n",
    "                slc.noiseCali.append(match[0])\n",
    "                print('Found noise calibration files: ', slc.noiseCali)\n",
    "\n",
    "def sort_caliFiles(slc):\n",
    "    if len(slc._tiffSrc)>0:\n",
    "        radioCaliList = []\n",
    "        noiseCaliList = []\n",
    "        s1_pat = re.compile(\".*(S1[AB]_.+?_\\d{4}\\d{2}\\d{2}T.*).zip.*\")\n",
    "        for swath in range(len(slc._tiffSrc)):\n",
    "            match = s1_pat.search(slc._tiffSrc[swath])\n",
    "            if match:\n",
    "                id_prefix=match.groups(1)[0]\n",
    "            else:\n",
    "                raise Exception('Unable to extract date in multi-slice scene: %s' %  slc._tiffSrc[swath])\n",
    "            for slc_ind in range(len(slc.radioCali)):\n",
    "                if id_prefix in slc.radioCali[slc_ind]:\n",
    "                    radioCaliList.append(slc.radioCali[slc_ind])\n",
    "       \t        if id_prefix in slc.noiseCali[slc_ind]:\n",
    "                    noiseCaliList.append(slc.noiseCali[slc_ind])\n",
    "        slc.radioCali=radioCaliList\n",
    "        slc.noiseCali=noiseCaliList        \n",
    "       \n",
    "    \n",
    "def write2flt( data, lat, lon, outName, nanvalue=-9999.90039062 ):\n",
    "    width = data.shape[1]\n",
    "    length = data.shape[0]\n",
    "    data[ (data == nanvalue) ] = 0\n",
    "    outm = np.matrix(data,np.float32)\n",
    "    gpmFile = open(outName, \"wb\")\n",
    "    outm.tofile(gpmFile)\n",
    "    gpmFile.close();\n",
    "\n",
    "def write2xml( data, xFirst, yFirst, dx, dy, outName, projection='lat/lon'):\n",
    "    length  = data.shape[0] \n",
    "    width   = data.shape[1] \n",
    "    xStep   = dx \n",
    "    yStep   = dy \n",
    "    proj    = projection\n",
    "    xmlName = outName + '.xml'\n",
    "    xmldict  = {'METADATA_LOCATION':xmlName,\n",
    "                'data_type':'Float',\n",
    "                'image_type':'BIL',\n",
    "                'Coordinate1':{'size':width,'startingValue':xFirst,'delta':xStep},\n",
    "                'Coordinate2':{'size':length,'startingValue':yFirst,'delta':yStep},\n",
    "                'FILE_NAME':outName,\n",
    "                'number_bands':1    }\n",
    "\n",
    "    demImage = createImage()\n",
    "    demImage.init(xmldict)\n",
    "    demImage.renderHdr()\n",
    "\n",
    "def readXMLDict(FileName):\n",
    "    xml = XmlParser\n",
    "    xml = XmlParser.XmlParser()\n",
    "    tmpdict  = xml.parse(FileName)[0]\n",
    "    xmldict = dict((k.replace(\" \", \"\").lower(), value) for (k, value) in tmpdict.items())\n",
    "    if xmldict == \"\":\n",
    "        print(\"No dictionay loaded from \"+FileName+\". \\n\");\n",
    "    else:\n",
    "        return(xmldict);    \n",
    "\n",
    "    \n",
    "'''\n",
    "To run radiometric and thermal noise calibration on a zip file, using vh-pol data, within a given extent:\n",
    "\n",
    "read_calibration_slc.py -zip S1A_IW_SLC__1SDV_20150315T231319_20150315T231349_005049_006569_0664.zip -od 20150315 -ext 34.6 34.65 -79.08 -78.97 -o -p vh -t noise\n",
    "\n",
    "To run radiometric calibration only using the topsApp.xml and scene xml file (using default vv-pol):\n",
    "\n",
    "read_calibration_slc.py -i topsApp.py -is s1_20150315.xml -o -t radio \n",
    "\n",
    "To output the radiometric calibration file only\n",
    "\n",
    "read_calibration_slc.py -zip S1A_IW_SLC__1SDV_20150315T231319_20150315T231349_005049_006569_0664.zip -od 20150315 -oc\n",
    "\n",
    "To output the radiometric and thermal noise calibration file (radio + noise) within a given extent\n",
    "\n",
    "read_calibration_slc.py -zip S1A_IW_SLC__1SDV_20150315T231319_20150315T231349_005049_006569_0664.zip -od 20150315 -o -ext 34.6 34.65 -79.08 -78.97 -t noise\n",
    "\n",
    "'''\n",
    "#read_calibration_slc.py( -zip {} -ext {} -od {} -o -t noise -n '{}'\\n\".format(ref_zip, bbox, reference_dir, ref_swath)\n",
    "def read_calibration_slc(fin, ext, odir, swath_num, type=\"noise\", output=True, pol=\"vv\", ck=False, oc=False): \n",
    "\n",
    "    if not isinstance(ext, list):\n",
    "        ext = ext.split()\n",
    "    ext = [np.float(ext[x]) for x in range(len(ext))]    \n",
    "    if type not in ('radio','noise'):\n",
    "        error('ERROR: type (-t) needs to be either radio or noise')\n",
    "    try:\n",
    "        swathList = swath_num.split()\n",
    "    except:\n",
    "        swathList = [1, 2, 3]\n",
    "\n",
    "    for swath in swathList:\n",
    "        slc = TOPS.createSentinel1()\n",
    "        slc.polarization = pol\n",
    "        slc.safe = fin\n",
    "        slc.swathNumber = int(swath)\n",
    "        #slc.product = TOPS.createTOPSSwathSLCProduct()\n",
    "        slc.regionOfInterest = ext\n",
    "        slc.product.bursts = iscesys.Component.createTraitSeq('bursts')\n",
    "        slc.output = os.path.join(odir,'IW{0}'.format(swath))\n",
    "        try:\n",
    "            slc.parse()\n",
    "        except Exception as err:\n",
    "            print('Could not extract swath {0} from {1}'.format(swath, slc.safe))\n",
    "            print('Generated: ', err)\n",
    "            continue\n",
    "        #slc.parse()\n",
    "        if len(slc.product.bursts) == 0:\n",
    "            continue\n",
    "        slc.radioCali = []\n",
    "        slc.noiseCali = []\n",
    "        locateCaliFile(slc,'radio',pol) \n",
    "        locateCaliFile(slc,'noise',pol) \n",
    "        sort_caliFiles(slc)\n",
    "        if slc._numSlices == 1:\n",
    "            slc.radioCali = (slc.product.numberOfBursts) * [slc.radioCali[0]]\n",
    "            slc.noiseCali = (slc.product.numberOfBursts) * [slc.noiseCali[0]]\n",
    "    \n",
    "        ################# Output calibrated SLC filesi (-o option)  ######################\n",
    "        if output:  # generate raw and calibrated slc\n",
    "\n",
    "            width  = slc._burstWidth\n",
    "            length = slc._burstLength\n",
    "  \n",
    "            ####Check if aux file corrections are needed\n",
    "            useAuxCorrections = False\n",
    "            if ('002.36' in slc.product.processingSoftwareVersion) and (slc.auxFile is not None):\n",
    "                useAuxCorrections = True\n",
    "  \n",
    "            ###If not specified, for single slice, use width and length from first burst\n",
    "            if width is None:\n",
    "                width = slc.product.bursts[0].numberOfSamples\n",
    "  \n",
    "            if length is None:\n",
    "                length = slc.product.bursts[0].numberOfLines\n",
    "\n",
    "            if os.path.isdir(slc.output):\n",
    "                print('Output directory exists. Overwriting ...')\n",
    "            else:\n",
    "                print('Creating directory {0} '.format(slc.output))\n",
    "                os.makedirs(slc.output)\n",
    "\n",
    "            prevTiff  = None\n",
    "            prevRadio = None\n",
    "            prevNoise = None\n",
    "            for index, burst in enumerate(slc.product.bursts):\n",
    "\n",
    "                print('############################################')\n",
    "                print('#  burst {} of IW{} '.format(index+1,swath))\n",
    "                print('############################################')\n",
    "                ####tiff for single slice\n",
    "                if (len(slc._tiffSrc) == 0) and (len(slc.tiff)==1):\n",
    "                    tiffToRead  = slc.tiff[0]\n",
    "                    radioToRead = slc.radioCali[0]\n",
    "                    noiseToRead = slc.noiseCali[0] \n",
    "                else: ##tiffSrc for multi slice\n",
    "                    tiffToRead  = slc._tiffSrc[index]\n",
    "                    radioToRead = slc.radioCali[index]\n",
    "                    noiseToRead = slc.noiseCali[index]\n",
    "                ###To minimize reads and speed up \n",
    "                if tiffToRead != prevTiff:\n",
    "                    src=None\n",
    "                    band=None\n",
    "                    src = gdal.Open(tiffToRead, gdal.GA_ReadOnly)\n",
    "                    fullWidth = src.RasterXSize\n",
    "                    fullLength = src.RasterYSize\n",
    "                    band = src.GetRasterBand(1)\n",
    "                    prevTiff  = tiffToRead\n",
    "                    xlist  = list(range(fullWidth))\n",
    "                    ylist  = list(range(fullLength))\n",
    "                    fullX, fullY = np.meshgrid(xlist, ylist)\n",
    "\n",
    "                #### Thermal Noise Calibration (Optional)\n",
    "                if ( noiseToRead != prevNoise ) and ( type == 'noise' ):\n",
    "                    if noiseToRead.startswith('/vsizip'): # read from zip file\n",
    "                       parts = noiseToRead.split(os.path.sep)\n",
    "                       if parts[2] == '':\n",
    "                           parts[2] = os.path.sep\n",
    "                       zipname = os.path.join(*(parts[2:-4]))\n",
    "                       fname = os.path.join(*(parts[-4:]))\n",
    "                       if not os.path.isfile(zipname):\n",
    "                           print('File ',zipname,' does not exist.')\n",
    "                           sys.exit(1)\n",
    "                       zf = zipfile.ZipFile(zipname, 'r')\n",
    "                       xmlstr = zf.read(fname)\n",
    "                    else: # Read out calibration from xml file\n",
    "                        fid = open(noiseToRead)\n",
    "                        xmlstr = fid.read()\n",
    "                    ### Extract the calibration LUT\n",
    "                    root   = ET.fromstring(xmlstr)\n",
    "                    try:\n",
    "                        caliRoot = root.find('noiseVectorList')\n",
    "                        nVectors  = int(caliRoot.items()[0][1])\n",
    "                    except:\n",
    "                        caliRoot = root.find('noiseRangeVectorList')\n",
    "                        nVectors  = int(caliRoot.items()[0][1])\n",
    "                    xxn = [];\n",
    "                    yyn = [];\n",
    "                    zzn = [];\n",
    "                    for child in caliRoot.getchildren():\n",
    "                        line  = int(child.find('line').text)\n",
    "                        pixel = list(map(int, child.find('pixel').text.split()))\n",
    "                        nPixel = int(child.find('pixel').items()[0][1])\n",
    "                        try:\n",
    "                            noiseLut = list(map(float, child.find('noiseLut').text.split()))\n",
    "                        except:\n",
    "                            noiseLut = list(map(float, child.find('noiseRangeLut').text.split()))\n",
    "                        xxn = xxn + pixel\n",
    "                        yyn = yyn + [line]*nPixel\n",
    "                        zzn = zzn + noiseLut\n",
    "                    if ck:\n",
    "                        if not any(zzn):  # if all-zeros in zzn\n",
    "                            zf.close()\n",
    "                            for izip in range(len(fin)):\n",
    "                                print(\"All-zero noise LUT in: {}\".format(fin[izip]))\n",
    "                                if(os.path.islink(fin[izip])):\n",
    "                                    print(\"{} is a symlink, doing unlink fo SLC with bad LUT\".format(fin[izip]))\n",
    "                                    os.unlink(fin[izip])\n",
    "                                else:\n",
    "                                    print('{} is not a symlink. Not removing file. Please run this again with symlinked files'.format(fin[izip]))\n",
    "                        else:\n",
    "                            print('LUT is okay for :', fin, '. Doing nothing.')\n",
    "                            sys.exit(1)\n",
    "                    else:    \n",
    "                        npt = len(zzn)\n",
    "                        coordn = np.hstack((np.array(xxn).reshape(npt,1),np.array(yyn).reshape(npt,1)))\n",
    "                        noise  = np.array(zzn).reshape(npt,1)\n",
    "                        print('Start 2D interpolation on noiseLut. This may take a while....')\n",
    "                        interpfn2 = interpnd(coordn,noise)\n",
    "                        noiseIntrp = interpfn2(fullX, fullY)\n",
    "                        prevNoise = noiseToRead\n",
    "                        if oc:\n",
    "                            ocfile = os.path.join(odir,'IW{0}'.format(swath),'noise-iw'+str(swath)+'.cal')\n",
    "                            write2flt( noiseIntrp, ylist, xlist, ocfile )\n",
    "                            write2xml( noiseIntrp, ylist[0], xlist[0], 1, 1, ocfile ) \n",
    "                            print('Writing noise calibration file to ', ocfile)\n",
    "\n",
    "                #### Radio Calibration\n",
    "                if radioToRead != prevRadio:\n",
    "                    ### Locate radiometric calibration xml file\n",
    "                    if radioToRead.startswith('/vsizip'): # read from zip file\n",
    "                        parts = radioToRead.split(os.path.sep)\n",
    "                        if parts[2] == '':\n",
    "                            parts[2] = os.path.sep\n",
    "                        zipname = os.path.join(*(parts[2:-4]))\n",
    "                        fname = os.path.join(*(parts[-4:]))\n",
    "                        if not os.path.isfile(zipname):\n",
    "                            print('File ',zipname,' does not exist.')\n",
    "                            sys.exit(1)\n",
    "                        zf = zipfile.ZipFile(zipname, 'r')\n",
    "                        xmlstr = zf.read(fname)\n",
    "                    else: # Read out calibration from xml file\n",
    "                        fid = open(radioToRead)\n",
    "                        xmlstr = fid.read()\n",
    "                    ### Extract the calibration LUT\n",
    "                    root   = ET.fromstring(xmlstr)\n",
    "                    caliRoot = root.find('calibrationVectorList')\n",
    "                    nVectors  = int(caliRoot.items()[0][1])\n",
    "                    xx = [];\n",
    "                    yy = [];\n",
    "                    zz = [];\n",
    "                    for child in caliRoot.getchildren():\n",
    "                        line  = int(child.find('line').text)\n",
    "                        pixel = list(map(int, child.find('pixel').text.split()))\n",
    "                        nPixel = int(child.find('pixel').items()[0][1])\n",
    "                        sigmaNought = list(map(float, child.find('sigmaNought').text.split()))\n",
    "                        xx = xx + pixel\n",
    "                        yy = yy + [line]*nPixel\n",
    "                        zz = zz + sigmaNought\n",
    "                    if ck:\n",
    "                        if not any(zz):  # if all-zeros in zz\n",
    "                            zf.close()\n",
    "                            for izip in range(len(fin)):\n",
    "                                os.unlink(fin[izip]) \n",
    "                            print('All-zero radiometric calibration LUT. Remove the symbolic link file ',fin )\n",
    "                    else: \n",
    "                        npt = len(zz)\n",
    "                        coord = np.hstack((np.array(xx).reshape(npt,1),np.array(yy).reshape(npt,1)))\n",
    "                        sigma  = np.array(zz).reshape(npt,1)\n",
    "                        print('Start 2D interpolation on sigmaNought. This may take a while....')\n",
    "                        interpfn1  = interpnd(coord,sigma)\n",
    "                        sigmaIntrp = interpfn1(fullX, fullY)\n",
    "                        prevRadio  = radioToRead\n",
    "                        if oc:   # read in calibration file\n",
    "                            ocfile = os.path.join(odir,'IW{0}'.format(swath),'radio-iw'+str(swath)+'.cal')\n",
    "                            write2flt( sigmaIntrp, ylist, xlist, ocfile )\n",
    "                            write2xml( sigmaIntrp, ylist[0], xlist[0], 1, 1, ocfile ) \n",
    "                            print('Writing radiometric calibration file to ', ocfile)\n",
    "\n",
    "                if ck:\n",
    "                    sys.exit(0)\n",
    "\n",
    "                #outfile1 = os.path.join(slc.output, 'burst_%02d_bk'%(index+1) + '.slc')\n",
    "                outfile2 = os.path.join(slc.output, 'burst_%02d'%(index+1) + '.slc')\n",
    " \n",
    "                ####Use burstnumber to look into tiff file\n",
    "                ####burstNumber still refers to original burst in slice\n",
    "                lineOffset = (burst.burstNumber-1) * burst.numberOfLines\n",
    "                data = band.ReadAsArray(0, lineOffset, burst.numberOfSamples, burst.numberOfLines)\n",
    "                datacrop = data[burst.firstValidLine:burst.lastValidLine, burst.firstValidSample:burst.lastValidSample]\n",
    " \n",
    "                ######Saving the radar cross-section\n",
    "                #if not os.path.isfile(outfile1):\n",
    "                #    fid = open(outfile1, 'wb')\n",
    "                #    outdata  = np.zeros((length,width), dtype=np.complex64)  #output radar cross-section\n",
    "                #    sigmaraw = np.square(abs(datacrop))\n",
    "                #    outdata[burst.firstValidLine:burst.lastValidLine, burst.firstValidSample:burst.lastValidSample] = sigmaraw\n",
    "                #    outdata.tofile(fid)\n",
    "                #    fid.close()\n",
    "                #    ####Render ISCE XML\n",
    "                #    slcImage = isceobj.createSlcImage()\n",
    "                #    slcImage.setByteOrder('l')\n",
    "                #    slcImage.setFilename(outfile1)\n",
    "                #    slcImage.setAccessMode('read')\n",
    "                #    slcImage.setWidth(width)\n",
    "                #    slcImage.setLength(length)\n",
    "                #    slcImage.setXmin(0)\n",
    "                #    slcImage.setXmax(width)\n",
    "                #    slcImage.renderHdr()\n",
    "                #    burst.image = slcImage\n",
    "                #else:\n",
    "                #    print('File {0} already exist. Skip file output.'.format(outfile1))\n",
    "    \n",
    "                fid = open(outfile2, 'wb')\n",
    "                outdata = np.zeros((length,width), dtype=np.complex64)  # convert to radar cross-section domain\n",
    "                calicrop = np.squeeze(sigmaIntrp[burst.firstValidLine:burst.lastValidLine, burst.firstValidSample:burst.lastValidSample])\n",
    "                sigmacali = np.square(abs(datacrop.real + 1j*datacrop.imag)/calicrop)\n",
    "                if type == 'radio':\n",
    "                    recomb = np.sqrt(sigmacali)*np.exp(1j*np.angle(datacrop))\n",
    "                elif type == 'noise':  #This means radio + noise\n",
    "                    noiseCrop = np.squeeze(noiseIntrp[burst.firstValidLine:burst.lastValidLine, burst.firstValidSample:burst.lastValidSample])\n",
    "                    noiseCorr = noiseCrop/np.square(calicrop)\n",
    "                    [ind0,ind1]   = np.where( abs(datacrop) == 0 )\n",
    "                    sigmadn       = sigmacali - noiseCorr\n",
    "                    sigmadn[ind0,ind1] = 0\n",
    "                    [ind2,ind3]   = np.where( sigmadn < 0 )\n",
    "                    sigmadn[ind2,ind3] = 1e-9  #instead of clipping at 0, assign a minimal value\n",
    "                    recomb    = np.sqrt(abs(sigmadn))*np.exp(1j*np.angle(datacrop))\n",
    "                outdata[burst.firstValidLine:burst.lastValidLine, burst.firstValidSample:burst.lastValidSample] = recomb\n",
    "                #Skip the correction for the Elevation Antenna Pattern for now \n",
    "                outdata.tofile(fid)\n",
    "                fid.close()\n",
    "                ####Render ISCE XML\n",
    "                slcImage = isceobj.createSlcImage()\n",
    "                slcImage.setByteOrder('l')\n",
    "                slcImage.setFilename(outfile2)\n",
    "                slcImage.setAccessMode('read')\n",
    "                slcImage.setWidth(width)\n",
    "                slcImage.setLength(length)\n",
    "                slcImage.setXmin(0)\n",
    "                slcImage.setXmax(width)\n",
    "                slcImage.renderHdr()\n",
    "            src=None\n",
    "            band=None\n",
    "\n",
    "def run_stackSlcDn_run2_5():\n",
    "    run1File = \"./run_files/run_01_unpack_topo_reference\"\n",
    "    run2File = \"./run_files/run_02_unpack_secondary_slc\"\n",
    "\n",
    "\n",
    "\n",
    "    reference_dir = os.path.join(wd, \"reference\")\n",
    "    secondary_dir = os.path.join(wd, \"secondary\")\n",
    "    refConfig = get_config_file(run1File, \"reference\")\n",
    "    print(refConfig)\n",
    "    \n",
    "    ref_zip = get_config_val(refConfig, 'dirname')\n",
    "    print(ref_zip)\n",
    "    ref_swath = get_config_val(refConfig, 'swaths')\n",
    "    print(ref_swath)\n",
    "    \n",
    "\n",
    "    try:\n",
    "        #read_calibration_slc(fin, ext, odir, swath_num, type=\"noise\", output=True, pol=\"vv\", ck=False, oc=False): \n",
    "        read_calibration_slc(ref_zip, bbox, reference_dir, ref_swath)\n",
    "        with open(run2File, 'r') as fp:\n",
    "            line = fp.readline()\n",
    "            while line:\n",
    "                if 'secondary' in line:\n",
    "                    secConfig = line.split(' ')[-1].strip()\n",
    "                    sec_date = get_date_from_str(os.path.basename(secConfig))\n",
    "                    print(\"{} : {}\".format(sec_date, os.path.basename(secConfig)))\n",
    "                    sec_zip = get_config_val(secConfig, 'dirname')\n",
    "                    sec_swath = get_config_val(secConfig, 'swaths')\n",
    "                    odir=os.path.join(secondary_dir, sec_date )\n",
    "                    print(\"{} : {} :{}\".format(sec_zip, sec_swath, odir ))\n",
    "                    read_calibration_slc(sec_zip, bbox, odir, sec_swath)\n",
    "            \n",
    "                line = fp.readline()\n",
    "    except Exception as e:\n",
    "        print(str(e))\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "MINLAT, MAXLAT, MINLON, MAXLON, MINLAT_LO, MAXLAT_HI, MINLON_LO, MAXLON_HI =get_bbox(ctx)\n",
    "print(\"{} {} {} {} {} {} {} {}\".format(MINLAT, MAXLAT, MINLON, MAXLON, MINLAT_LO, MAXLAT_HI, MINLON_LO, MAXLON_HI))\n",
    "\n",
    "bbox = \"{} {} {} {}\".format(MINLAT, MAXLAT, MINLON, MAXLON)\n",
    "print(bbox)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Download Sentinel-1 data SLC**\n",
    "\n",
    "The SLCs identified in the parameters cell are downloaded. Note that these files are significant in size (around 5GB each) which will collectively require a significant amount of time to download (expect around 2 minutes per SLC)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "path = SLC_PATH\n",
    "download_slcs(path)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "***Download Orbit Files based on SLC***\n",
    "\n",
    "The orbit files identified based on the set of SLCs downloaded in the previous step. They are placed in the orbits subdirectory."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "get_orbit_files()\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Download DEM data and generate DEM file**\n",
    "Similarly to the orbit files, the DEM information to retain is computed based on the bounding box.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "download_dem()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Fix paths in dem file**\n",
    "\n",
    "This step invokes the isce2 application fixImageXML.py in order to fill in the complete path to the files referenced in the dem file.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "import os\n",
    "cwd = os.getcwd()\n",
    "print(\"cwd : {}\".format(cwd))\n",
    "if os.path.exists(\"dem.txt\"):\n",
    "    cmd = [\"awk\", \"'/wgs84/ {print $NF;exit}'\", \"dem.txt\"]\n",
    "    WGS84 = run_cmd_output(cmd).decode(\"utf-8\").strip()\n",
    "    wgs84_file = os.path.join(cwd, WGS84)\n",
    "    print(\"WGS84 : a{}b\".format(wgs84_file))\n",
    "    if os.path.exists(wgs84_file):\n",
    "        print(\"Found wgs84 file: {}\".format(wgs84_file))\n",
    "        fix_cmd = [\"{}/applications/fixImageXml.py\".format(ISCE_HOME), \"--full\", \"-i\", \"{}\".format(wgs84_file) ]\n",
    "        run_cmd(fix_cmd) \n",
    "    else:\n",
    "        print(\"NO WGS84 FILE FOUND : {}\".format(wgs84_file))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### AUX_CAL file download ####\n",
    "\n",
    "The following calibration auxliary (AUX_CAL) file is used for **antenna pattern correction** to compensate the range phase offset of SAFE products with **IPF verison 002.36** (mainly for images acquired before March 2015). If all your SAFE products are from another IPF version, then no AUX files are needed. Check [ESA document](https://earth.esa.int/documents/247904/1653440/Sentinel-1-IPF_EAP_Phase_correction) for details. \n",
    "\n",
    "Make sure your aws credentials are fresh (i.e. run 'aws-login -p default' in a terminal window) and run the command below to download the AUX_CAL file once and store it somewhere (_i.e._ ~/aux/aux_cal) so that you can use it all the time, for `stackSentinel.py -a` or `auxiliary data directory` in `topsApp.py`.\n",
    "\n",
    "```aws s3 cp --recursive s3://nisar-dev-ondemand/S1_aux/  ~/aux/aux_cal/```\n",
    "\n",
    "In the cell below, the necessary AUX_CAL file is copied to the AuxDir directory from this download."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%bash\n",
    "mkdir -p ./AuxDir\n",
    "cp ~/aux/aux_cal/S1A/S1A_AUX_CAL_V20140915T100000_G20151125T103928.SAFE/data/s1a-aux-cal.xml ./AuxDir"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "The scripts provides support for Sentinel-1 TOPS stack processing. Currently supported workflows include a coregistered stack of SLC, interferograms, offsets, and coherence. \n",
    "\n",
    "`stackSentinel.py` generates all configuration and run files required to be executed on a stack of Sentinel-1 TOPS data. When stackSentinel.py is executed for a given workflow (-W option) a **configs** and **run_files** folder is generated. No processing is performed at this stage. Within the run_files folder different run\\_#\\_description files are contained which are to be executed as shell scripts in the run number order. Each of these run scripts call specific configure files contained in the “configs” folder which call ISCE in a modular fashion. The configure and run files will change depending on the selected workflow. To make run_# files executable, change the file permission accordingly (e.g., `chmod +x run_01_unpack_slc`).\n",
    "\n",
    "```bash\n",
    "stackSentinel.py -H     #To see workflow examples,\n",
    "stackSentinel.py -h     #To get an overview of all the configurable parameters\n",
    "```\n",
    "\n",
    "Required parameters of stackSentinel.py include:\n",
    "\n",
    "```cfg\n",
    "-s SLC_DIRNAME          #A folder with downloaded Sentinel-1 SLC’s. \n",
    "-o ORBIT_DIRNAME        #A folder containing the Sentinel-1 orbits. Missing orbit files will be downloaded automatically\n",
    "-a AUX_DIRNAME          #A folder containing the Sentinel-1 Auxiliary files\n",
    "-d DEM_FILENAME         #A DEM (Digital Elevation Model) referenced to wgs84\n",
    "```\n",
    "\n",
    "In the following, different workflow examples are provided. Note that stackSentinel.py only generates the run and configure files. To perform the actual processing, the user will need to execute each run file in their numbered order.\n",
    "\n",
    "In all workflows, coregistration (-C option) can be done using only geometry (set option = geometry) or with geometry plus refined azimuth offsets through NESD (set option = NESD) approach, the latter being the default. For the NESD coregistrstion the user can control the ESD coherence threshold (-e option) and the number of overlap interferograms (-O) to be used in NESD estimation.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%bash\n",
    "rm -rf ./run_files "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "master_date=get_master_date(ctx)\n",
    "print(\"master_date : {}\".format(master_date))\n",
    "\n",
    "if master_date:\n",
    "    print(\"MASTER_DATE exists:\".format(master_date) )   \n",
    "    cmd = [\n",
    "        \"/opt/isce2/src/isce2/contrib/stack/topsStack/stackSentinel.py\",  \"-s\", \"zip/\", \"-d\", \"{}\".format(wgs84_file), \"-a\", \"AuxDir/\", \"-m\", \"{}\".format(master_date), \"-o\", \"orbits\", \n",
    "        \"-b\", \"\\\"{} {} {} {}\\\"\".format(MINLAT, MAXLAT, MINLON, MAXLON), \n",
    "        \"-W\", \"slc\", \"-C\", \"geometry\"\n",
    "    ]\n",
    "else:\n",
    "    print(\"MASTER_DATE DOES NOT EXIST\")          \n",
    "    cmd = [\n",
    "        \"/opt/isce2/src/isce2/contrib/stack/topsStack/stackSentinel.py\",  \"-s\", \"zip/\", \"-d\", \"{}\".format(wgs84_file), \"-a\", \"AuxDir/\", \"-o\", \"orbits\", \n",
    "        \"-b\", \"\\\"{} {} {} {}\\\"\".format(MINLAT, MAXLAT, MINLON, MAXLON), \n",
    "        \"-W\", \"slc\", \"-C\", \"geometry\"\n",
    "    ]\n",
    "run_cmd(cmd)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**stackSlcDn_run2_5**\n",
    "\n",
    "An intermediate run file must be manually added between run files 2 and 3 in order to perform radiometric and thermal noise cancellation. This step generates this run file, executing against each SLC file"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "stackSlcDn_run2_5()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**run_01_unpack_slc_topo_reference:**\n",
    "Unpacks the reference SLC files using ISCE readers, also unpackaging the antenna elevation pattern correction file if necessary. This run file also produces the reference geometry files that are consumed downstream.\n",
    "\n",
    "Includes commands to unpack Sentinel-1 TOPS SLCs using ISCE readers. For older SLCs which need antenna elevation pattern correction, the file is extracted and written to disk. For newer version of SLCs which don’t need the elevation antenna pattern correction, only a gdal virtual “vrt” file (and isce xml file) is generated. The “.vrt” file points to the Sentinel SLC file and reads them whenever required during the processing. If a user wants to write the “.vrt” SLC file to disk, it can be done easily using gdal_translate (e.g. gdal_translate –of ENVI File.vrt File.slc). \n",
    "The “run_01_unpack_slc_topo_reference” also includes a command that refers to the config file of the stack reference, which includes configuration for running topo for the stack reference. Note that in the pair-wise processing strategy one should run topo (mapping from range-Doppler to geo coordinate) for all pairs. However, with stackSentinel, topo needs to be run only one time for the reference in the stack. \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%bash\n",
    "echo \"## STEP 1 ##\"\n",
    "\n",
    "start=`date +%s`\n",
    "echo \"sh run_files/run_01_unpack_topo_reference\"\n",
    "sh run_files/run_01_unpack_topo_reference\n",
    "end=`date +%s`\n",
    "runtime1=$((end-start))\n",
    "echo \"STEP 1 RUN TIME : $runtime1\""
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**run_02_unpack_secondary_slc:**\n",
    "\n",
    "Unpack Secondary SLCs\n",
    "In a manner similar to the SLCs in Run File 01 above, this run file unpacks the secondary SLCs from each of the input SLC zip files."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%bash\n",
    "echo \"## STEP 2 ##\"\n",
    "\n",
    "start=`date +%s`\n",
    "Num=`cat run_files/run_02_unpack_secondary_slc | wc | awk '{print $1}'`\n",
    "echo $Num\n",
    "echo \"cat run_files/run_02_unpack_secondary_slc | parallel -j2 --eta --load 50%\"\n",
    "cat run_files/run_02_unpack_secondary_slc | parallel -j2 --eta --load 50%\n",
    "end=`date +%s`\n",
    "\n",
    "runtime2=$((end-start))\n",
    "echo STEP 2 RUN TIME : $runtime2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### run_02.5_slc_noise_calibration ###\n",
    "This run file, which was manually introduced in step 19, runs radiometric and thermal noise calibration using vh-pol data against the SLCs."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%bash\n",
    "echo \"## STEP 2.5 ##\"\n",
    "\n",
    "start=`date +%s`\n",
    "echo \"run_files/run_02.5_slc_noise_calibration| parallel -j2 --eta --load 50%\"\n",
    "cat run_files/run_02.5_slc_noise_calibration|parallel -j2 --eta --load 50%\n",
    "end=`date +%s`\n",
    "\n",
    "runtime2x5=$((end-start))\n",
    "echo STEP 2.5 RUN TIME : $runtime2x5\n",
    "    "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**run_03_average_baseline:**\n",
    "\n",
    "Computes average baseline for the stack. These baselines are not used for processing anywhere. They are only an approximation and can be used for plotting purposes."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%bash\n",
    "echo \"## STEP 3 ##\"\n",
    "\n",
    "start=`date +%s`\n",
    "echo \"cat run_files/run_03_average_baseline | parallel -j2 --eta --load 50%\"\n",
    "cat run_files/run_03_average_baseline | parallel -j2 --eta --load 50%\n",
    "end=`date +%s`\n",
    "runtime3=$((end-start))\n",
    "echo STEP 3 RUN TIME : $runtime3"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**STEP 4 : run_04_fullBurst_geo2rdr:**\n",
    "\n",
    "Full Burst geo2rdr\n",
    "This run file estimates geometrical offsets between secondary burst overlaps and the stack reference burst overlaps. The secondaries are then resampled to the stack reference burst overlaps.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%bash\n",
    "echo \"## STEP 4 ##\"\n",
    "\n",
    "start=`date +%s`\n",
    "\n",
    "echo \"cat run_files/run_04_fullBurst_geo2rdr  | parallel -j2 --eta --load 50%\"\n",
    "cat run_files/run_04_fullBurst_geo2rdr  | parallel -j2 --eta --load 50%\n",
    "end=`date +%s`\n",
    "runtime4=$((end-start))\n",
    "echo STEP 4 RUN TIME : $runtime4"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**STEP 5 : run_05_fullBurst_resample:**\n",
    "\n",
    "Using orbit and DEM data, this run file computes geometrical offsets among all secondary SLCs and the stack reference. These offsets, with the misregistration time series are used for precise coregistration of each burst SLC."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%bash\n",
    "echo \"## STEP 5 ##\"\n",
    "\n",
    "start=`date +%s`\n",
    "echo \"cat run_files/run_05_fullBurst_resample  | parallel -j2 --eta --load 50%\"\n",
    "cat run_files/run_05_fullBurst_resample  | parallel -j2 --eta --load 50%\n",
    "end=`date +%s`\n",
    "runtime5=$((end-start))\n",
    "echo STEP 5 RUN TIME : $runtime5\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**STEP 6 : run_06_timeseries_misreg:**\n",
    "\n",
    "This run file extracts the valid region between burst SLCs at the overlap area of the bursts. This region changes slightly for different acquisitions and must be retained for merging the bursts to eliminate invalid data.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%bash\n",
    "echo \"## STEP 6 ##\"\n",
    "\n",
    "start=`date +%s`\n",
    "echo \"sh run_files/run_06_extract_stack_valid_region\"\n",
    "sh run_files/run_06_extract_stack_valid_region\n",
    "end=`date +%s`\n",
    "runtime6=$((end-start))\n",
    "echo STEP 6 RUN TIME : $runtime6"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**run_07_geo2rdr_resample: Merge Reference Secondary SLCs**\n",
    "\n",
    "This run file merges all bursts for the reference and coregistered SLCs. Geometry files are also merged.\n",
    "Using orbit and DEM, geometrical offsets among all secondary SLCs and the stack reference is computed. The geometrical offsets, together with the misregistration time-series (from previous step) are used for precise coregistration of each burst SLC. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%bash\n",
    "echo \"## STEP 7 ##\"\n",
    "\n",
    "start=`date +%s`\n",
    "\n",
    "FILE=run_files/run_07_merge\n",
    "if test -f \"$FILE\"; then\n",
    "    echo \"cat run_files/run_07_merge  | parallel -j2 --eta --load 50%\"\n",
    "    cat run_files/run_07_merge | parallel -j2 --eta --load 50%\n",
    "else\n",
    "    echo \"cat run_files/run_07_merge_reference_secondary_slc  | parallel -j2 --eta --load 50%\"\n",
    "    cat run_files/run_07_merge_reference_secondary_slc | parallel -j2 --eta --load 50%\n",
    "fi\n",
    "\n",
    "end=`date +%s`\n",
    "runtime7=$((end-start))\n",
    "echo STEP 7 RUN TIME : $runtime7"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Create Dataset**\n",
    "The final step ‘packages’ up the results into a single sub-directory named as coregistered_slcs-<start_date>-<end_date> and generates json files describing the dataset (e.g. the polygon  coordinates of the bounding box) and associated metadata (e.g. included SLCs)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "create_dataset(bbox)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "run_stop_time = get_current_time()\n",
    "print(\"PGE run stop time : {}\".format(run_stop_time))\n",
    "total_run_time = run_stop_time - run_start_time\n",
    "print(\"Total Run Time = {}\".format(total_run_time))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "celltoolbar": "Tags",
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.9.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}