{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 风云4A(FY4A) AGRI L1数据读取基础介绍"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1. 支持的文件类型"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "支持NetCDF格式的风云4A(FY4A) AGRI数据(全圆盘和中国区)\n",
    "\n",
    "示例:\n",
    "\n",
    "全圆盘:\n",
    "    \n",
    "    FY4A-_AGRI--_N_DISK_1047E_L1-_FDI-_MULT_NOM_20190807060000_20190807061459_4000M_V0001.HDF\n",
    "\n",
    "中国区:\n",
    "\n",
    "    FY4A-_AGRI--_N_REGC_1047E_L1-_FDI-_MULT_NOM_20190807045334_20190807045750_1000M_V0001.HDF\n",
    "\n",
    "*全圆盘标识：DISK , 中国区标识：REGC.*\n",
    "\n",
    "数据链接:\n",
    "\n",
    "    实时数据(30天内)及文档:\n",
    "        https://fy4.nsmc.org.cn/data/en/data/realtime.html\n",
    "    历史数据 (2018-03-12 -- ):\n",
    "        http://satellite.nsmc.org.cn/PortalSite/Data/Satellite.aspx\n",
    "    FY4A官方应用平台:\n",
    "        http://rsapp.nsmc.org.cn/geofy/"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2. 定标"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "有以下三种:\n",
    "\n",
    "1. 原始数字量化值 (所有通道)\n",
    "\n",
    "2. 反射率 (C01 - C06)\n",
    "\n",
    "3. 辐射率和亮温 (C07 - C14)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3. 例子"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 安装"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "```\n",
    "$ conda install -c conda-forge satpy\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 加载数据"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import os, glob\n",
    "from satpy.scene import Scene"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "['C01',\n",
       " 'C02',\n",
       " 'C03',\n",
       " 'C04',\n",
       " 'C05',\n",
       " 'C06',\n",
       " 'C07',\n",
       " 'C08',\n",
       " 'C09',\n",
       " 'C10',\n",
       " 'C11',\n",
       " 'C12',\n",
       " 'C13',\n",
       " 'C14',\n",
       " 'satellite_azimuth_angle',\n",
       " 'satellite_zenith_angle',\n",
       " 'solar_azimuth_angle',\n",
       " 'solar_glint_angle',\n",
       " 'solar_zenith_angle']"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# 加载FY4A文件\n",
    "filenames = glob.glob('/xin/data/FY4A/20190807/FY4A-_AGRI*4000M_V0001.HDF')\n",
    "\n",
    "# 创建scene对象\n",
    "scn = Scene(filenames, reader='agri_l1')\n",
    "\n",
    "# 查看可用的通道\n",
    "scn.available_dataset_names()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 以红外通道为例\n",
    "ir_channel = 'C12'\n",
    "scn.load([ir_channel], generate=False, calibration='brightness_temperature')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 全圆盘图"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 在notebook中显示\n",
    "scn.show(ir_channel)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<img src='./figures/agri_C12.png'>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 保存到文件\n",
    "# scn.save_dataset(ir_channel, filename='{sensor}_{name}.png')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 真彩色全圆盘图"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Too many possible datasets to load for DatasetID(name=None, wavelength=3.9, resolution=None, polarization=None, calibration=None, level=None, modifiers=())\n",
      "Too many possible datasets to load for 3.9\n",
      "Too many possible datasets to load for 3.9\n",
      "Too many possible datasets to load for 3.9\n",
      "Too many possible datasets to load for DatasetID(name=None, wavelength=3.9, resolution=None, polarization=None, calibration=None, level=None, modifiers=())\n",
      "Too many possible datasets to load for 3.9\n",
      "Too many possible datasets to load for 3.9\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "['ash',\n",
       " 'dust',\n",
       " 'fog',\n",
       " 'green',\n",
       " 'green_snow',\n",
       " 'ir108_3d',\n",
       " 'ir_cloud_day',\n",
       " 'natural_color',\n",
       " 'natural_color_sun',\n",
       " 'night_background',\n",
       " 'night_background_hires',\n",
       " 'overview',\n",
       " 'overview_sun',\n",
       " 'true_color']"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# 查看可用的合成选项\n",
    "scn.available_composite_names()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Required file type 'agri_l1_4000m_geo' not found or loaded for 'solar_zenith_angle'\n",
      "Required file type 'agri_l1_4000m_geo' not found or loaded for 'satellite_azimuth_angle'\n",
      "Required file type 'agri_l1_4000m_geo' not found or loaded for 'satellite_zenith_angle'\n",
      "Required file type 'agri_l1_4000m_geo' not found or loaded for 'solar_azimuth_angle'\n"
     ]
    }
   ],
   "source": [
    "# 注：这步需要大内存 (取决于cpu核数)\n",
    "# 可以查看FAQ中关于内存的讨论：\n",
    "#    https://satpy.readthedocs.io/en/latest/faq.html\n",
    "\n",
    "composite = 'true_color'\n",
    "scn.load([composite])\n",
    "scn.show(composite)\n",
    "# scn.save_dataset(composite, filename='{sensor}_{name}.png')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<img src='./figures/agri_true_color.png'>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 特定区域图"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "以下以台风利奇马为例。\n",
    "\n",
    "首先，我们需要定义地图投影和区域，然后将数据投影到该区域上。\n",
    "\n",
    "我们可以用`Pyresample`来定义区域，并将其写入到`area.yaml`方便以后直接调用。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "from pyresample import get_area_def\n",
    "\n",
    "area_id = 'lekima'\n",
    "\n",
    "x_size = 549\n",
    "y_size = 499\n",
    "area_extent = (-1098006.560556, -967317.140452, 1098006.560556, 1026777.426728)\n",
    "projection = '+proj=laea +lat_0=19.0 +lon_0=128.0 +ellps=WGS84'\n",
    "description = \"Typhoon Lekima\"\n",
    "proj_id = 'laea_128.0_19.0'\n",
    "\n",
    "areadef = get_area_def(area_id, description, proj_id, projection,x_size, y_size, area_extent)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "你可以用[coord2area_def.py](https://github.com/pytroll/satpy/blob/main/utils/coord2area_def.py)程序来直接生成区域定义。\n",
    "\n",
    "比如用`python coord2area_def.py lekima_4km laea 10 28 118 138 4`可得到之前定义的利奇马区域：\n",
    "\n",
    "```\n",
    "lekima_4km:\n",
    "  description: lekima_4km\n",
    "  projection:\n",
    "    proj: laea\n",
    "    ellps: WGS84\n",
    "    lat_0: 19.0\n",
    "    lon_0: 128.0\n",
    "  shape:\n",
    "    height: 499\n",
    "    width: 549\n",
    "  area_extent:\n",
    "    lower_left_xy: [-1098006.560556, -967317.140452]\n",
    "    upper_right_xy: [1098006.560556, 1026777.426728]\n",
    "```\n",
    "\n",
    "然后将该定义拷贝到`$PPP_CONFIG_DIR/areas.yaml`中，即可直接调用。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 如果你已经添加区域到areas.yaml中，可直接调用:\n",
    "os.environ['PPP_CONFIG_DIR'] = '/yin_raid/xin/satpy_config/'\n",
    "lekima_scene = scn.resample('lekima_4km')\n",
    "\n",
    "# 否则需要使用之前定义的区域:\n",
    "# lekima_scene = scn.resample(areadef)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "lekima_scene.show(composite)\n",
    "# lekima_scene.save_dataset(composite, filename='{sensor}_{name}_resampled.png')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<img src='./figures/agri_true_color_resampled.png'>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "如果想利用自定义的colormap来生成图像（如下图），请参阅关于`enhancement`的notebook。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<img src='./figures/agri_C12_resampled_colorize.png'>"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "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.7.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}