1batch_topsApp.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# batch running topsApp"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Setting basic config and parameters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#########################################################\n",
    "#  Modify following variables according to your condition\n",
    "# ========================================================\n",
    "# Folder contains insar data\n",
    "sar_dir = '/data/sentinel1/data'\n",
    "# DEM path \n",
    "dem = '/data/sentinel1/DEM/demLat_N38_N40_Lon_E097_E101.dem'\n",
    "# dem = None # auto download\n",
    "# Folder contains Precise orbits\n",
    "orbit = '/data/sentinel1/aux/aux_poeorb'\n",
    "# Folder contains Auxiliary data\n",
    "aux = '/data/sentinel1/aux/aux_cal'\n",
    "\n",
    "range_looks = 20\n",
    "azimuth_looks = 4\n",
    "filter_strength = 0.01\n",
    "swath_number = None\n",
    "unwrapper_name = 'snaphu'\n",
    "\n",
    "roi = [38.4573, 39.0911, 98.6675, 99.3582]\n",
    "# [S, N, W, E] convention (lat/lon)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The code below are functions about generating ifg pairs and topsApp.xml file"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import csv\n",
    "import logging\n",
    "import logging.config\n",
    "import os\n",
    "from datetime import datetime\n",
    "from pathlib import Path\n",
    "from tqdm import tqdm\n",
    "from isceobj.XmlUtil import FastXML as xml\n",
    "\n",
    "logging.config.fileConfig(os.path.join(os.getcwd(), 'logging.conf'))\n",
    "logger = logging.getLogger('batch_tops')\n",
    "\n",
    "\n",
    "def generate_topsApp_xml(\n",
    "        reference_dir, secondary_dir, reference_file, secondary_file, roi, xml_file):\n",
    "    '''导出干涉对xml配置文件'''\n",
    "\n",
    "    # reference and secondary information\n",
    "    fix_info = {\"orbit directory\": orbit,\n",
    "                \"auxiliary data directory\": aux}\n",
    "    if swath_number is not None:\n",
    "        fix_info.update({\"swath number\": fix_info})\n",
    "\n",
    "    reference = {\"safe\": reference_file,\n",
    "                 \"output directory\": reference_dir}\n",
    "    reference.update(fix_info)\n",
    "\n",
    "    secondary = {\"safe\": secondary_file,\n",
    "                 \"output directory\": secondary_dir}\n",
    "    secondary.update(fix_info)\n",
    "\n",
    "    tops = xml.Component('topsinsar')\n",
    "    # parameters\n",
    "    tops[\"Sensor name\"] = 'SENTINEL1'\n",
    "\n",
    "    tops['reference'] = reference\n",
    "    tops['secondary'] = secondary\n",
    "\n",
    "    if roi:\n",
    "        tops[\"region of interest\"] = roi\n",
    "    if dem:\n",
    "        tops[\"demFilename\"] = Path(dem).name\n",
    "    tops[\"use GPU\"] = True\n",
    "    tops[\"range looks\"] = range_looks\n",
    "    tops[\"azimuth looks\"] = azimuth_looks\n",
    "    tops[\"do unwrap\"] = True\n",
    "    tops[\"unwrapper name\"] = unwrapper_name\n",
    "    tops[\"filter strength\"] = filter_strength\n",
    "\n",
    "    tops.writeXML(xml_file, root='topsApp')\n",
    "\n",
    "\n",
    "def mapper_date_data(folder, suffix):\n",
    "    '''映射日期到数据文件\n",
    "    Select the data files from the folder,and then return \n",
    "    a mapper from year to data file\n",
    "    '''\n",
    "    files_s1 = sorted([i for i in os.listdir(folder)\n",
    "                      if i.split('.')[-1] == suffix])\n",
    "    dates = [i[17:25] for i in files_s1]\n",
    "    mapper = dict()\n",
    "    for i, date in enumerate(dates):\n",
    "        mapper.update({date: files_s1[i]})\n",
    "    return mapper\n",
    "\n",
    "\n",
    "def load_interferogram_pairs(path_list, delimiter=','):\n",
    "    '''加载干涉对文件'''\n",
    "    with open(path_list) as f:\n",
    "        ms_list = list(csv.reader(f, delimiter=delimiter))\n",
    "    return ms_list\n",
    "\n",
    "\n",
    "def day_interval(date_start, date_end):\n",
    "    date_start = datetime.strptime(date_start, '%Y%m%d')\n",
    "    date_end = datetime.strptime(date_end, '%Y%m%d')\n",
    "    interval = (date_end - date_start).days\n",
    "\n",
    "    return interval\n",
    "\n",
    "\n",
    "def generate_interferogram_pairs(sar_dir,  pair_file, sar_suffix='zip',\n",
    "                                 max_interval=2, max_day=180, delimiter=','):\n",
    "    '''导出干涉对文件。\n",
    "    Generate and save a reference-secondary pair automatically. Time interval \n",
    "    of reference-secondary pair is 1 and 2\n",
    "\n",
    "    dates: str \n",
    "        date string with format of '%Y%m%d'\n",
    "    max_interval: int\n",
    "        干涉对最大相邻获取间隔, 相邻间隔为1\n",
    "    max_day:int\n",
    "        干涉对最大的间隔天数\n",
    "    pair_file: str\n",
    "        the path of file out.if None, will not save. \n",
    "    delimiter: str\n",
    "        reference与secondary日期之间的分割符号\n",
    "    '''\n",
    "    mapper = mapper_date_data(sar_dir, suffix=sar_suffix)\n",
    "    dates = sorted(mapper.keys())\n",
    "\n",
    "    num = len(dates)\n",
    "    reference_secondary_pair_list = []\n",
    "    for i, date in enumerate(dates):\n",
    "        dti_temp = 1\n",
    "        while dti_temp <= max_interval:\n",
    "            if i + dti_temp < num:\n",
    "                if day_interval(date, dates[i + dti_temp]) < max_day:\n",
    "                    reference_secondary_pair_list.append(\n",
    "                        [date, dates[i + dti_temp]])\n",
    "                    dti_temp += 1\n",
    "                else:\n",
    "                    break\n",
    "            else:\n",
    "                break\n",
    "\n",
    "    with open(pair_file, 'w', newline='') as f:\n",
    "        csv_writer = csv.writer(f, delimiter=delimiter)\n",
    "        csv_writer.writerows(reference_secondary_pair_list)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Firstly, let's generate interferogram pairs and save those pairs to a file. You can modify the ifg pairs manually after generation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "################ generate interferogram pairs ################\n",
    "\n",
    "# the folder to store interferograms\n",
    "home_dir = Path('/data/sentinel1/result')\n",
    "# the file contains ifg pairs\n",
    "pair_file = home_dir / 'ifg_pairs.csv'\n",
    "\n",
    "generate_interferogram_pairs(sar_dir, pair_file, sar_suffix='zip',\n",
    "                             max_interval=3, max_day=180, delimiter=',')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now, let's load designed ifg pairs and batch running topsApp.py"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# %%   batch run topsApp.py\n",
    "ifg_pairs = load_interferogram_pairs(pair_file)\n",
    "mapper = mapper_date_data(sar_dir, suffix='zip')\n",
    "\n",
    "for i, ifg in enumerate(tqdm(ifg_pairs)):\n",
    "    reference, secondary = ifg\n",
    "    reference_file = os.path.join(sar_dir, mapper[reference])\n",
    "    secondary_file = os.path.join(sar_dir, mapper[secondary])\n",
    "\n",
    "    ifg_name = f'{reference}_{secondary}'\n",
    "    ifg_dir = home_dir / f'interferogams/{ifg_name}'\n",
    "    if not ifg_dir.is_dir():\n",
    "        ifg_dir.mkdir(parents=True)\n",
    "\n",
    "    # create xml file\n",
    "    reference_dir = ifg_dir / f'{reference}'\n",
    "    secondary_dir = ifg_dir / f'{secondary}'\n",
    "    xml_file = ifg_dir / 'topsApp.xml'\n",
    "\n",
    "    generate_topsApp_xml(reference_dir, secondary_dir,\n",
    "                         reference_file, secondary_file, roi, xml_file)\n",
    "\n",
    "    # cd ifg folder\n",
    "    os.chdir(ifg_dir)\n",
    "\n",
    "    # link dem file to current folder\n",
    "    print(ifg_name)\n",
    "    cmd_link = f'ln -sf {str(Path(dem).parent)}/* ./'\n",
    "    os.system(cmd_link)\n",
    "\n",
    "    # execute topsApp.py\n",
    "    cmd_tops = f'topsApp.py {xml_file} --start=startup'\n",
    "    logger.info(cmd_tops)\n",
    "    os.system(cmd_tops)\n"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3.7.13 ('geo')",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "name": "python",
   "version": "3.7.13"
  },
  "orig_nbformat": 4,
  "vscode": {
   "interpreter": {
    "hash": "630114ded1cb6136faf7350823e03790fe4102d781074abef124d6f86b9fe410"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}