Thule-mix-LAADS.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "0fd1f287",
   "metadata": {
    "tags": []
   },
   "source": [
    "# **MODIS vs THAAO comparison of Integrated Water Vapor (IWV)**\n",
    "[MODIS](https://modis.gsfc.nasa.gov/about/specifications.php): Moderate Resolution Imaging Spectroradiometer  <br>\n",
    "[THAAO](https://www.thuleatmos-it.it/index.php): Thule High Arctic Atmospheric Observatory, Greenland  (coordinates  __76.5145, -68.7432__, altitude a.s.l. __220 m__)  <br>\n",
    "\n",
    "__MXD05_L2 IWV product__ <br>\n",
    "_References_ <br>\n",
    "[MYD05_L2](https://atmosphere-imager.gsfc.nasa.gov/products/water-vapor) <br>\n",
    "\n",
    "_Algorithm Theoretical Basis Documents (ATBDs)_ <br>\n",
    "[ATBD](https://atmosphere-imager.gsfc.nasa.gov/documentation/atbds-plans-guides) <br>\n",
    "[The MODIS Near-IR Water Vapor Algorithm](https://atmosphere-imager.gsfc.nasa.gov/sites/default/files/ModAtmo/atbd_mod03_0.pdf) <br>\n",
    "[MODIS Atmospheric Profile Retrieval\n",
    "ATBD MOD07 Collection 6](https://modis-images.gsfc.nasa.gov/_docs/MOD07_atbd_v7_April2011.pdf) <br>\n",
    "[Atmospheric Correction Algorithm: Spectral\n",
    "Reflectances MOD08](https://modis.gsfc.nasa.gov/data/atbd/atbd_mod08.pdf) <br>\n",
    "[MODIS Atmosphere L3 Gridded Product\n",
    "(ATBD) & Users Guide](https://atmosphere-imager.gsfc.nasa.gov/sites/default/files/ModAtmo/L3_ATBD_C6_C61_2019_02_20.pdf) <br>\n",
    "\n",
    "__THAAO IWV product__ <br>\n",
    "_References_ <br>\n",
    "[The water Vapor Emission Spectrometer for Polar Atmospheres at 22 GHz](https://www.thuleatmos-it.it/data/vespa22/index.php)<br> \n",
    "[Integrated Water Vapor](https://www.thuleatmos-it.it/dataaccess/vapor/index.php)<br>\n",
    " "
   ]
  },
  {
   "cell_type": "raw",
   "id": "c807a143-1c94-4588-b8eb-9148faa2fc66",
   "metadata": {},
   "source": [
    "%run reader.ipynb"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8d72e112-5d1d-492b-be8a-776c8243a0b4",
   "metadata": {
    "tags": []
   },
   "source": [
    "## __Data filtering main parameters: IWV, distance from THAAO, deltatime__ \n",
    "Parameters used for MODIS data selection based on maximum reasonable IWV value (max_IWV), and maximum distance of the MODIS pixel nearest to the THAAO station. <br>\n",
    "The deltatime parameter (expressed in minutes) sets the maximum time difference between a MODIS  measure and a VESPA measure for a refined comparison. <br>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "66719500-7130-4c22-9185-85dd9cdff4d7",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# passaggio parametri di filtraggio dati\n",
    "%store -r filter_param\n",
    "\n",
    "max_IWV = filter_param.get('max_IWV')\n",
    "max_IWV_unit = filter_param.get('max_IWV_unit')\n",
    "max_dist_from_THAAO = filter_param.get('max_dist_from_THAAO')\n",
    "max_dist_from_THAAO_unit = filter_param.get('max_dist_from_THAAO_unit')\n",
    "deltatime = filter_param.get('deltatime')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "14886aaf-4d85-4144-82e2-647b389cfe1e",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# print di controllo\n",
    "print(f'max_IWV = {max_IWV} mm\\n')\n",
    "print(f'max_dist_from_THAAO = {max_dist_from_THAAO} Km\\n')\n",
    "print(f'deltatime = {deltatime}\\n')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "acda0539-af9d-4cd5-bf34-c9f5f283a389",
   "metadata": {},
   "source": [
    "## __Loading MXD05_L2 and THAAO IWV time series__ \n",
    "\n",
    "Satellite MODIS MXD05_L2 data products covering the Greenland region on a time interval spanning from July 2016 to January 2022 were first downloaded from the LAADS repository by means of custom Python scripts. Timeseries of IWV IR and ancillary data products were then extracted on a point geometry centered on the THAAO Observatory, Greenland from the MXD05_L2 products and packed together as csv files for further analisys. <br>\n",
    "\n",
    "Thule (Greenland) IWV data timeseries have been downloaded from the [Integrated Water Vapor](https://www.thuleatmos-it.it/dataaccess/vapor/index.php) THAAO data access facility were the data are freely accessible. Likewise MODIS satellite IWV product, these data contain vertical column integrated amount  water vapor or precipitable water retrieved by means of the water Vapor Emission Spectrometer for Polar Atmospheres at 22 GHz (VESPA-22).  <br> \n",
    "\n",
    "The available sensors' time series are readed by the PreProc.ipynb\n",
    "IWV timeseries of MODIS in .csv and THAAO-VESPA-22 in .txt formats are read into pandas dataframes. <br>\n",
    "MODIS Near-IR product has an estimated additional 13% error column derived from the maximum error claimed in the NearIR IWV product ATBD (see sect. 3.1.4, page 16). <br>\n",
    "THAAO IWV data are given in $kg/m^2$ units and thus divided for factor of 10 to have homogeoeous units with MODIS IWV data that are given in $cm$. <br>  "
   ]
  },
  {
   "cell_type": "raw",
   "id": "bc84c1f6-8fcb-42ce-955e-4f238e2f274b",
   "metadata": {},
   "source": [
    "mod_data = modis_read('MODIS/MOD05_5km_totale.csv')\n",
    "myd_data = modis_read('MODIS/MYD05_5km_totale.csv')\n",
    "vespa_data = vespa_read('VESPA/Vapor_20160716_20220311.txt')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ac7f8320",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "%store -r mod_data \n",
    "%store -r myd_data \n",
    "%store -r mxd_data \n",
    "%store -r vespa_data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "54ef9893-3c7b-4d04-af98-f1b4d12e5941",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "vespa_data.tail(2)     # ok: raw per test"
   ]
  },
  {
   "cell_type": "raw",
   "id": "14f28345-7696-4aa8-b70e-d33994abd991",
   "metadata": {},
   "source": [
    "display(vespa_data.tail(2))\n",
    "print('aaaaaaaaaaaaaaaaa')"
   ]
  },
  {
   "cell_type": "raw",
   "id": "d8fe6e1b-4f4f-4d1c-940f-8da0315e065e",
   "metadata": {
    "tags": []
   },
   "source": [
    "myd_data.tail(2)       # ok: raw per test"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "81e5c1cb-1a9c-44ee-9dd5-c0439d209cdf",
   "metadata": {},
   "source": [
    "#### **Controllo numero di punti mod, myd e mxd per max IWV e distanza da THAAO selezionati**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "51659f2b",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "number_mod_data_p = len(mod_data.IWV_MODIS) \n",
    "number_myd_data_p = len(myd_data.IWV_MODIS)\n",
    "number_mxd_data_p = len(mxd_data.IWV_MODIS)\n",
    "number_mxd_data_points = number_mod_data_p + number_myd_data_p"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f8fdac7b-81d4-4875-b3f1-e4a070da26ed",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "print(f'max IWV = {max_IWV} [mm]')\n",
    "print(f'max distance from THAOO = {max_dist_from_THAAO} [km]\\n')\n",
    "print(f'Total number of MOD data points selected =    {number_mod_data_p}')\n",
    "print(f'Total number of MYD data points selected =    {number_myd_data_p}\\n')\n",
    "print(f'Total number of MODIS data points selected = {number_mxd_data_p}\\n')\n",
    "print(f'Total number of MODIS data points selected = {number_mxd_data_points}\\n')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1e1f771e-13b1-4b92-ab8a-0013de2deb44",
   "metadata": {},
   "source": [
    "#### **Controllo distanza e lat,lon dei pixel MODIS dal THAAO**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "cb3e2bfb-33b5-4b03-bffc-12b32947bae4",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "print(f'MOD distance [km]:  min = {mod_data.diff_distance.min()}   max = {mod_data.diff_distance.max()}')\n",
    "print(f'MYD distance [km]:  min = {myd_data.diff_distance.min()}   max = {myd_data.diff_distance.max()}\\n')\n",
    "print(f'MOD latitude  :  min = {mod_data.latitude.min()}    max = {mod_data.latitude.max()}')\n",
    "print(f'MOD longitude :  min = {mod_data.longitude.min()}   max = {mod_data.longitude.max()} \\n')\n",
    "print(f'MYD latitude  :  min = {myd_data.latitude.min()}    max = {myd_data.latitude.max()}')\n",
    "print(f'MYD longitude :  min = {myd_data.longitude.min()}   max = {myd_data.longitude.max()} \\n')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f082439b-8ea2-421a-b075-80dc3e9cf025",
   "metadata": {},
   "source": [
    "#### **Numero punti per mese** "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "35744cc0-f19b-4693-bdeb-6db55df696b4",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "mod_data.Time.groupby([mod_data['Time'].dt.month.rename('month')]).agg({'count'}).transpose()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4ebb04fd-9a29-4700-a063-8a79610237d4",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "myd_data.Time.groupby([myd_data['Time'].dt.month.rename('month')]).agg({'count'}).transpose()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "0ede1b8a-efc7-4576-9ffc-df0f0bf4deae",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "mxd_data.Time.groupby([mxd_data['Time'].dt.month.rename('month')]).agg({'count'}).transpose()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6803ddf7-544a-4c19-a5dd-b01b55e31955",
   "metadata": {},
   "source": [
    "## **Other sensors**"
   ]
  },
  {
   "cell_type": "raw",
   "id": "7eb9de12-825a-444a-b8b7-365864866873",
   "metadata": {
    "tags": []
   },
   "source": [
    "\n",
    "wdir_hatpro = 'HATPRO/'\n",
    "DATA_dir = '/home/luca/WORK/DATA/Groundbased/THULE/'\n",
    "ELAB_dir = '/home/luca/WORK/ELAB/Groundbased/THULE/'\n",
    "\n",
    "\n",
    "def read_hatpro_year(fname):\n",
    "    year = fname[-15:-11]\n",
    "    return year\n",
    "import glob\n",
    "from pathlib import *\n",
    "\n",
    "p = Path('.')\n",
    "butta = list(p.glob('**/*.ipynb'))\n",
    "\n",
    "p_data = PurePath(p.home(),'WORK','DATA','Groundbased','THULE')\n",
    "p_elab = PurePath(p.home(),'WORK','ELAB','Groundbased','THULE')\n",
    "#year = [read_hatpro_year(y) for y in fnlist]\n",
    "butta\n",
    "#p_elab.mkdir(mode=0o777, parents=False)\n",
    "dir = p_elab"
   ]
  },
  {
   "cell_type": "raw",
   "id": "c3b66da4",
   "metadata": {
    "tags": []
   },
   "source": [
    "def set_IO_dirs(path_segments):\n",
    "    \"\"\" Final path_segments as a list, to be added to\n",
    "    the baseline directory fo a complete absolute path\n",
    "    to DATA and ELAB dirs.\n",
    "    \"\"\"\n",
    "    from pathlib import Path\n",
    "    \n",
    "    # writing the path segments\n",
    "    homedir = Path.home()          # similar to ~/\n",
    "    data_dir_segments = [homedir,'WORK','DATA','Groundbased','THULE']\n",
    "    [data_dir_segments.append(s) for s in path_segments]\n",
    "    #elab_dir_segments = data_dir_segments\n",
    "    #elab_dir_segments = data_dir_segments[:]\n",
    "    import copy\n",
    "    elab_dir_segments = copy.deepcopy(data_dir_segments)\n",
    "    \n",
    "    elab_dir_segments[2] = 'ELAB'\n",
    "    return data_dir_segments, elab_dir_segments\n",
    "\n",
    "def check_dirs_existence_mkdir(p_data, p_elab):\n",
    "    \"\"\" Check out DATA and ELAB I/O dirs existence and \n",
    "    eventually creates them.\n",
    "    \"\"\"\n",
    "    p_data.mkdir(parents=True, exist_ok=True)\n",
    "    p_elab.mkdir(parents=True, exist_ok=True)\n",
    "    return\n",
    "\n",
    "def field_campaign_name():\n",
    "    pass\n",
    "    return fc_name\n",
    "\n",
    "def sensor_name():\n",
    "    pass\n",
    "    return sensor_name\n",
    "\n",
    "\n",
    "\n",
    "from pathlib import *\n",
    "\n",
    "#transforming path segments into concrete paths\n",
    "path_segments  = ['HATPRO']\n",
    "path_data, path_elab = set_IO_dirs(path_segments)\n",
    "\n",
    "print(f'data: {path_data}    id: {id(path_data)}')\n",
    "print(f'elab: {path_elab}    id: {id(path_elab)}')\n",
    "\n",
    "#p_data, p_elab = Path(*final_path_segments)\n",
    "p_data = Path(*path_data)\n",
    "p_elab = Path(*path_elab)\n",
    "\n",
    "print(f'data: {p_data}    id: {id(p_data)}')\n",
    "print(f'elab: {p_elab}    id: {id(p_elab)}')\n",
    "\n",
    "check_dirs_existence_mkdir(p_data, p_elab)"
   ]
  },
  {
   "cell_type": "raw",
   "id": "2fe523e3",
   "metadata": {
    "tags": []
   },
   "source": [
    "year=[2016,2017,2018,2019]\n",
    "#gnss_data = pd.DataFrame()\n",
    "hatpro_data = pd.DataFrame()\n",
    "\n",
    "for year in year:\n",
    "#    gnss_data = gnss_data.append(gnss_read(\n",
    "#        f'GNSS/THU3_{year}global.txt',year))\n",
    "    hatpro_data = hatpro_data.append(hatpro_read(\n",
    "        f'HATPRO/IWV_{year}_HATPRO.txt',year))"
   ]
  },
  {
   "cell_type": "raw",
   "id": "15ed6000",
   "metadata": {
    "tags": []
   },
   "source": [
    "# HATPRO = hatpro_read('/home/luca/WORK/DATA/Groundbased/THULE/HATPRO/IWV_2016_HATPRO.txt')\n",
    "\n",
    "#gnss_data.head()\n",
    "hatpro_data.head()"
   ]
  },
  {
   "cell_type": "raw",
   "id": "be0f40d1",
   "metadata": {
    "tags": []
   },
   "source": [
    "# old directories\n",
    "HATPRO = hatpro_read('/home/luca/WORK/DATA/Groundbased/THULE/HATPRO/IWV_2017_HATPRO.txt')\n",
    "VESPA = vespa_read('/home/luca/WORK/DATA/Groundbased/THULE/VESPA22/Vapor_20160716_20220311.txt')\n",
    "MOD_IR = modis_read('/home/luca/WORK/DATA/MODIS/THULE/Thule_LAADS_IWV_IR/MOD05_5km_totale.csv')\n",
    "MYD_IR = modis_read('/home/luca/WORK/DATA/MODIS/THULE/Thule_LAADS_IWV_IR/MYD05_5km_totale.csv')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0492d15c",
   "metadata": {
    "tags": []
   },
   "source": [
    "## **MODIS and THAAO Integrated Water Vapor (IWV) total timeseries**\n",
    "Simple linear plots of the whole MODIS and THAAO IVW timeseries."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f3ab002d-eb0d-4744-93e9-12c095cb9a43",
   "metadata": {
    "tags": []
   },
   "source": [
    "#### **Visualization options**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "53c403ae-470a-4a6f-b05c-70ada27d461b",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "width = 1100\n",
    "height = 300 \n",
    "size = 3\n",
    "scale = 0.3\n",
    "color_MOD = 'cyan'\n",
    "color_MYD = '#F2583E'\n",
    "color_THAAO = 'green'\n",
    "\n",
    "alpha=100\n",
    "dynspread=False\n",
    "datashade=False\n",
    "legend='top_right'\n",
    "\n",
    "# text options\n",
    "xlabel='Time'\n",
    "ylabel='IWV    [mm]'\n",
    "\n",
    "plot_defaults = dict(xlabel=xlabel, ylabel=ylabel,height=height,\n",
    "                width=width,alpha=alpha,dynspread=dynspread,\n",
    "                datashade=datashade,size=size, legend=legend)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c313ff8d-1298-4add-bac3-d9a5a02a3282",
   "metadata": {
    "tags": [],
    "toc-hr-collapsed": true
   },
   "source": [
    "### **Time-series with MODIS data selected by max IWV and distance from THAAO**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "979a8198-6b73-43f0-8cb4-4069d921a3af",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# TODO capire perché bisogna importare di nuovo hvplot\n",
    "import hvplot.pandas"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c5276e8c-3bdf-4dee-91ac-0ffd4f5b3e50",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "mod_timeseries = mod_data.IWV_MODIS.hvplot.scatter(color=color_MOD, label='MOD',\n",
    "                                                 title='IWV MOD',\n",
    "                                                 **plot_defaults)\n",
    "myd_timeseries = myd_data.IWV_MODIS.hvplot.scatter(color=color_MYD, label='MYD',\n",
    "                                                 title='IWV MYD',\n",
    "                                                 **plot_defaults)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6f165c40-e83b-4493-a8fe-58754b541f58",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "thaao_timeseries = vespa_data.IWV_THAAO.hvplot.scatter(color=color_THAAO, label='VESPA',\n",
    "                                                       fields={'date': 'Date'},\n",
    "                                                       title='IWV THAAO - VESPA', \n",
    "                                                       **plot_defaults)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3ad7d5da-9937-45c1-9a9d-a8baf03b3430",
   "metadata": {},
   "source": [
    "### **Time-series plots with errors**"
   ]
  },
  {
   "cell_type": "raw",
   "id": "149776d8-eb39-4516-b035-911b1af5ea1a",
   "metadata": {
    "tags": []
   },
   "source": [
    "# MODIS data with errors     # ok: raw per test\n",
    "(mod_timeseries*mod_data.hvplot.errorbars(y='IWV_MODIS', yerr1='Err_IWV_MODIS', hover=False)\n",
    "+ myd_timeseries*myd_data.hvplot.errorbars(y='IWV_MODIS', yerr1='Err_IWV_MODIS', hover=False)\n",
    "+ (mod_timeseries*mod_data.hvplot.errorbars(y='IWV_MODIS', yerr1='Err_IWV_MODIS', hover=False)\n",
    "   *myd_timeseries*myd_data.hvplot.errorbars(y='IWV_MODIS', yerr1='Err_IWV_MODIS', hover=False)\n",
    "  ).opts(title='MXD')).cols(1)"
   ]
  },
  {
   "cell_type": "raw",
   "id": "5b4d4e00-539e-424c-b7bf-4e9fbc941d86",
   "metadata": {
    "tags": []
   },
   "source": [
    "# VESPA data with errors     # ok: raw per test\n",
    "print(f'VESPA   min_err = {vespa_data.err_IWV.min()},  max_err = {vespa_data.err_IWV.max()}')\n",
    "thaao_timeseries*vespa_data.hvplot.errorbars(y='IWV_THAAO', yerr1='err_IWV', hover=False)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "183b6aef-6e0d-4efc-89f3-a9bba88638b6",
   "metadata": {},
   "source": [
    "### **Time-series plots**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1727babb-67c4-4411-92b6-ff5f70d38797",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "((mod_timeseries*myd_timeseries).opts(title='MXD') \n",
    " + thaao_timeseries).cols(1)"
   ]
  },
  {
   "cell_type": "raw",
   "id": "5fe31d5a-8a35-41cc-8db5-34216c15e734",
   "metadata": {
    "tags": []
   },
   "source": [
    "mod_data['Err_IWV_MODIS'].head(1)     # ok: raw per test"
   ]
  },
  {
   "cell_type": "raw",
   "id": "cc596a31-101f-4703-9a79-2f2f09563338",
   "metadata": {},
   "source": [
    "vespa_data.head(1)     # ok: raw per test"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4c7c6e1f",
   "metadata": {
    "tags": [],
    "toc-hr-collapsed": true
   },
   "source": [
    "### **Time-series plots of all MODIS and VESPA-22 IWV data on the same time reference** <br>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "affdb0b9-c87e-4faf-93d2-03d208d0f3e4",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "title_mod=f'IWV MOD - VESPA {chr(0x20)*50} # MOD data = {number_mod_data_p}'\n",
    "title_myd=f'IWV MYD - VESPA {chr(0x20)*50} # MYD data = {number_myd_data_p}'\n",
    "title_mxd=f'IWV MXD - VESPA {chr(0x20)*50} # MXD data = {number_mxd_data_points}'"
   ]
  },
  {
   "cell_type": "raw",
   "id": "f5a5a08b-07a5-45c6-9e58-c069c317a9af",
   "metadata": {
    "tags": []
   },
   "source": [
    "# MODIS MOD and MYD time series      # ok: raw per test\n",
    "((mod_timeseries *thaao_timeseries).opts(title=title_mod) \n",
    "+ (myd_timeseries * thaao_timeseries).opts(title=title_myd)).cols(1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "609993cf",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "(mod_timeseries*myd_timeseries * thaao_timeseries).opts(title=title_mxd, show_grid=True)"
   ]
  },
  {
   "cell_type": "raw",
   "id": "06366dda-9d10-4a90-a096-bbd6fd4f2b40",
   "metadata": {
    "tags": []
   },
   "source": [
    "# MODIS and VESPA data with MODIS errors     # ok: raw per test\n",
    "\n",
    "(mod_timeseries*myd_timeseries * thaao_timeseries\n",
    "* mod_data.hvplot.errorbars(y='IWV_MODIS', yerr1='Err_IWV_MODIS', hover=False)\n",
    "* myd_data.hvplot.errorbars(y='IWV_MODIS', yerr1='Err_IWV_MODIS', hover=False)\n",
    ").opts(title=title_mxd)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ee428428",
   "metadata": {
    "tags": []
   },
   "source": [
    "### __QA__ \n",
    "### __TODO: inserire il discorso sul check dei Quality Assessment: tutti i dati sono best quality__ <br>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "27150287-b092-40a9-810b-f7148d25d4e8",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# processing and visualization parameters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "faf9599b-6dda-485c-9cb6-db37240e20fb",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "%store -r new_mod_vespa \n",
    "%store -r new_myd_vespa \n",
    "%store -r new_mxd_vespa "
   ]
  },
  {
   "cell_type": "raw",
   "id": "2942fe90-4165-4a7d-9193-c2d44a56a587",
   "metadata": {},
   "source": [
    "new_myd_vespa = closest_in_time(myd_data,vespa_data,'IWV_MODIS','IWV_THAAO',on='Time',delta_time=deltatime)\n",
    "new_mod_vespa = closest_in_time(mod_data,vespa_data,'IWV_MODIS','IWV_THAAO',on='Time',delta_time=deltatime)"
   ]
  },
  {
   "cell_type": "raw",
   "id": "7dce41e7-7df2-48ba-90a9-252a09d3e4b5",
   "metadata": {
    "tags": []
   },
   "source": [
    "new_mxd_vespa = pd.concat([new_mod_vespa,new_myd_vespa], sort=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c60fd43c-c932-46f0-b210-7a3e4dd0fd77",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "new_mod_vespa.head(1)     # ok: raw per test"
   ]
  },
  {
   "cell_type": "raw",
   "id": "73339421-6e0f-4f87-947d-d8fb149e2cdd",
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    },
    "tags": []
   },
   "source": [
    "new_myd_vespa.tail(1)     # ok: raw per test"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "82f6a692-abb7-4cf7-ac3d-cf9963930e7f",
   "metadata": {},
   "source": [
    "#### **Controllo numero di punti mod, myd e mxd per il deltatime selezionato**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8ea118e3-9fb5-4523-a26b-e330bd4aa28d",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "num_mod_data_p_Dt = len(new_mod_vespa.IWV_MODIS) \n",
    "num_myd_data_p_Dt = len(new_myd_vespa.IWV_MODIS)\n",
    "num_mxd_data_p_Dt = len(new_mxd_vespa.IWV_MODIS)\n",
    "num_modmyd_data_p_Dt = num_mod_data_p_Dt + num_myd_data_p_Dt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "bed71abd-d68e-4abb-9ad7-d4ff058aacf0",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "print(f'Number of data points selected with deltatime = {deltatime}in')\n",
    "print(f'within a max distance from THAOO = {max_dist_from_THAAO} km\\n')\n",
    "print(f'Selected number of MOD data points =          {num_mod_data_p_Dt}')\n",
    "print(f'Selected number of MYD data points =          {num_myd_data_p_Dt}\\n')\n",
    "print(f'Total selected number of MXD data points     =  {num_mxd_data_p_Dt}')\n",
    "print(f'Total selected number of MOD+MYD data points =  {num_modmyd_data_p_Dt}\\n')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1e8091c8",
   "metadata": {},
   "source": [
    "#### **Controllo distanza e lat,lon dei pixel MODIS dal THAAO**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "49d3788e",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "print(f'MOD distance [km]:  min = {new_mod_vespa.diff_distance.min()}   max = {new_mod_vespa.diff_distance.max()}')\n",
    "print(f'MYD distance [km]:  min = {new_myd_vespa.diff_distance.min()}   max = {new_myd_vespa.diff_distance.max()}\\n')\n",
    "print(f'MOD latitude  :  min = {new_mod_vespa.latitude.min()}    max = {new_mod_vespa.latitude.max()}')\n",
    "print(f'MOD longitude :  min = {new_mod_vespa.longitude.min()}   max = {new_mod_vespa.longitude.max()} \\n')\n",
    "print(f'MYD latitude  :  min = {new_myd_vespa.latitude.min()}    max = {new_myd_vespa.latitude.max()}')\n",
    "print(f'MYD longitude :  min = {new_myd_vespa.longitude.min()}   max = {new_myd_vespa.longitude.max()} \\n')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "18cab918",
   "metadata": {},
   "source": [
    "#### **Numero punti per mese** "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a2871942",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "new_mod_vespa.Time.groupby([new_mod_vespa['Time'].dt.month.rename('month')]).agg({'count'}).transpose()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "847593d1",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "new_myd_vespa.Time.groupby([new_myd_vespa['Time'].dt.month.rename('month')]).agg({'count'}).transpose()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d468320f-157a-416c-977d-a8b12b3736f4",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "new_mxd_vespa.Time.groupby([new_mxd_vespa['Time'].dt.month.rename('month')]).agg({'count'}).transpose()"
   ]
  },
  {
   "cell_type": "raw",
   "id": "bb2a0d15",
   "metadata": {
    "tags": []
   },
   "source": [
    "newdf['RMSE_IWV']=((newdf.IWV_MODIS-newdf.IWV_THAAO)**2).mean()**.5\n",
    "newdf.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "86a6e44a-136b-4281-8d40-3bf8cc513689",
   "metadata": {},
   "source": [
    "#### **Controllo dei pixel selezionati con deltatime, max and min time difference tra le misure MODIS e VESPA**\n",
    "TODO"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "35da9a7e-a72b-4288-8463-d529ed253a69",
   "metadata": {},
   "source": [
    "#### **Controllo VZA mod e myd per il deltatime selezionato**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "fe76082b-f41e-4323-a57f-a7beb0d47a3d",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "print(f'MOD View Zenith Angle:  min = {new_mod_vespa.vza.min()}   max = {new_mod_vespa.vza.max()}')\n",
    "print(f'MYD View Zenith Angle:  min = {new_myd_vespa.vza.min()}   max = {new_myd_vespa.vza.max()}\\n')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "bed29998-1992-43b3-bdd2-0b3e5893fb2e",
   "metadata": {},
   "source": [
    "### **New time-series definition after time range selection**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "22537f15-cde6-414f-8a2b-7ff595742fce",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "t_mod_vespa=f'IWV MOD - VESPA{chr(0x20)*50}deltatime={deltatime}{chr(0x20)*50}# MOD data = {num_mod_data_p_Dt}'\n",
    "t_myd_vespa=f'IWV MYD - VESPA{chr(0x20)*50}deltatime={deltatime}{chr(0x20)*50}# MYD data = {num_myd_data_p_Dt}'\n",
    "t_mxd_vespa=f'IWV MXD - VESPA{chr(0x20)*50}deltatime={deltatime}{chr(0x20)*50}# MXD data = {num_mxd_data_p_Dt}'\n",
    "\n",
    "new_mod_vespa_ts = new_mod_vespa.IWV_MODIS.hvplot.scatter(x='Time_IWV_MODIS',\n",
    "                                                          color=color_MOD, \n",
    "                                                          label='MOD',\n",
    "                                                          title=t_mod_vespa,\n",
    "                                                          **plot_defaults)\n",
    "new_myd_vespa_ts = new_myd_vespa.IWV_MODIS.hvplot.scatter(x='Time_IWV_MODIS',\n",
    "                                                          color=color_MYD, \n",
    "                                                          label='MYD',\n",
    "                                                          title=t_myd_vespa,\n",
    "                                                          **plot_defaults)\n",
    "\n",
    "new_vespa_mod_ts = new_mod_vespa.IWV_THAAO.hvplot.scatter(x='Time_IWV_THAAO',\n",
    "                                                         color=color_THAAO,\n",
    "                                                         label='VESPA',                                                        \n",
    "                                                         **plot_defaults)\n",
    "new_vespa_myd_ts = new_myd_vespa.IWV_THAAO.hvplot.scatter(x='Time_IWV_THAAO',\n",
    "                                                         color=color_THAAO,\n",
    "                                                         label='VESPA',                                                        \n",
    "                                                         **plot_defaults)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c0630145",
   "metadata": {},
   "source": [
    "### **New time-series plots after time range selection** "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9b7c4fcd-5bd9-4bbb-8898-ef1e5cedc647",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "(new_mod_vespa_ts * new_myd_vespa_ts\n",
    "* new_vespa_mod_ts * new_vespa_myd_ts).opts(title=t_mxd_vespa, show_grid=True)"
   ]
  },
  {
   "cell_type": "raw",
   "id": "955571da-7033-4bcd-9d0f-95cbf7e10557",
   "metadata": {
    "tags": []
   },
   "source": [
    "# MODIS and VESPA data with MODIS errors     # ok: raw per test\n",
    "\n",
    "(new_mod_vespa_ts *new_myd_vespa_ts\n",
    "* new_mod_vespa.hvplot.errorbars(y='IWV_MODIS', yerr1='Err_IWV_MODIS', hover=False)\n",
    "* new_myd_vespa.hvplot.errorbars(y='IWV_MODIS', yerr1='Err_IWV_MODIS', hover=False)   \n",
    "* new_vespa_mod_ts * new_myd_vespa_ts).opts(title=t_mxd_vespa, show_grid=True)"
   ]
  },
  {
   "cell_type": "raw",
   "id": "2be52f04-2623-4136-867c-e5ae21d3bbde",
   "metadata": {
    "tags": []
   },
   "source": [
    "new_mxd_vespa.head(2)     # ok: raw per test"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "79c6907f-3a06-4121-9cf7-f4dfa7b8cb8f",
   "metadata": {},
   "source": [
    "# **Data Comparison MODIS vs VESPA-22**"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d4861cb5",
   "metadata": {},
   "source": [
    "### Regressione con Holoviews"
   ]
  },
  {
   "cell_type": "raw",
   "id": "079fa72b",
   "metadata": {
    "tags": []
   },
   "source": [
    "def regression_holoviews(x,y,color='orange', overlay_values=True):\n",
    "    \"\"\"\n",
    "    Input: x values and y values\n",
    "    overlay values: overlay of coefficients on the plot\n",
    "\n",
    "    Return a list having:\n",
    "    \n",
    "    1- all regress element\n",
    "    2- slope value,\n",
    "    3- intercept value \"\"\"\n",
    "\n",
    "    coords = {'x':x,'y':y}\n",
    "    xxxx= hv.Scatter(coords)\n",
    "    regress=hv.Slope.from_scatter(xxxx).opts(color=color)\n",
    "    \n",
    "    if overlay_values == True:      \n",
    "        return [regress*hv.Text(max(x)-0.3,min(y)+0.1, \n",
    "                                f'y = {regress.slope} x + {regress.y_intercept}'),\n",
    "                regress.slope,\n",
    "                regress.y_intercept]\n",
    "    else:\n",
    "        return [regress,\n",
    "                regress.slope,\n",
    "                regress.y_intercept]"
   ]
  },
  {
   "cell_type": "raw",
   "id": "76f4d015",
   "metadata": {
    "tags": []
   },
   "source": [
    "regress_MOD = regression_holoviews(newdf_MOD.IWV_MODIS,newdf_MOD.IWV_THAAO)\n",
    "regress_MYD = regression_holoviews(newdf_MYD.IWV_MODIS,newdf_MYD.IWV_THAAO,color='black')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e210e8ef",
   "metadata": {},
   "source": [
    "### **Regressione con sklearn**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "2031cf74-24f6-463e-97f4-0c061ec59c10",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# TODO da capire come passare lib e funzioni senza far rigirare i nb \n",
    "# (cioè facendo girare i nb separatamente, non attraverso Thule_main.ipynb)\n",
    "%run reader.ipynb\n",
    "%run functions.ipynb"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "486422cc",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "regress_mod_sk = regression_sklearn(new_mod_vespa.IWV_MODIS,new_mod_vespa.IWV_THAAO)\n",
    "regress_mod_sk_text = regression_eq_text(new_mod_vespa.IWV_MODIS,\n",
    "                                         new_mod_vespa.IWV_THAAO, regress_mod_sk)\n",
    "regress_myd_sk = regression_sklearn(new_myd_vespa.IWV_MODIS,\n",
    "                                    new_myd_vespa.IWV_THAAO, color='black')\n",
    "regress_myd_sk_text = regression_eq_text(new_myd_vespa.IWV_MODIS,\n",
    "                                         new_myd_vespa.IWV_THAAO, regress_myd_sk)\n",
    "regress_mxd_vespa_sk = regression_sklearn(new_mxd_vespa.IWV_MODIS,new_mxd_vespa.IWV_THAAO, \n",
    "                                    color='cyan')\n",
    "regress_mxd_vespa_sk_text = regression_eq_text(new_mxd_vespa.IWV_MODIS,\n",
    "                                         new_mxd_vespa.IWV_THAAO, regress_mxd_vespa_sk)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5ef27bb4",
   "metadata": {},
   "source": [
    "#### **Bisettrice**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a15187dc",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "slope = hv.Slope(1, 0).opts(line_color='red')   # generic bisector"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "55ef1dce",
   "metadata": {
    "tags": []
   },
   "source": [
    "## **MODIS vs THAAO Integrated Water Vapor (IWV) timeseries scatter plots**"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6bd89661-9003-4bee-9440-a2fba5a8ce3e",
   "metadata": {
    "tags": []
   },
   "source": [
    "## **Scatter plots**"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "db9d6f8a-2091-42d0-823a-99397290a37e",
   "metadata": {},
   "source": [
    "#### **Visualization options: scatter plots**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d1863810-4a11-4e9a-9857-cdd5757f221b",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# line-dots colors                \n",
    "color_MOD = 'cyan'\n",
    "color_MYD = '#F2583E'\n",
    "color_MXD = 'magenta'\n",
    "color_THAAO = 'green'\n",
    "\n",
    "# general default settings\n",
    "x = ''\n",
    "y = ''\n",
    "xlabel = 'IWV MOD    [mm]' \n",
    "ylabel = 'IWV VESPA    [mm]'\n",
    "color = ''\n",
    "label = ''\n",
    "height = 500\n",
    "width = 600\n",
    "c = ''\n",
    "s = ''\n",
    "cmap = ''\n",
    "alpha = 0.5\n",
    "scale = 1\n",
    "#size = 2\n",
    "legend = 'top_right'\n",
    "dynspread = False\n",
    "datashade = False  # datashading not available, può essere tolto\n",
    "grid=True\n",
    "title = ''\n",
    " \n",
    "# ATTENZIONE nell'uso di **kwargs (come splot_defaults): \n",
    "# - il parametro opzionale s che viene passato come kwds  \n",
    "#   NON può essere definito se non vieno utilizzato:\n",
    "#   le celle non girano. \n",
    "#   Va definito dopo con splot_defaults.update(s='colonna_df'), \n",
    "#   SOLO quando si usa\n",
    "# - il parametro size blocca l'effetto del parametro s \n",
    "#   (i punti hanno dim bloccata a size che non varia con s)\n",
    "\n",
    "splot_defaults = dict(x=x, y=x,\n",
    "                      xlabel=xlabel, ylabel=ylabel,\n",
    "                      color=color, label=label,\n",
    "                      height=height, width=width,\n",
    "                      c=c,\n",
    "                      cmap=cmap, \n",
    "                      alpha=alpha,\n",
    "                      scale=scale, \n",
    "                      legend=legend,\n",
    "                      dynspread=dynspread,datashade=datashade,\n",
    "                      grid=grid,\n",
    "                      title = title)\n",
    "\n",
    "## metodo \n",
    "#s_=''\n",
    "#splot_defaults.update(s=s)\n",
    "#print(f' metodo 2: {splot_defaults}')\n",
    "print(splot_defaults)"
   ]
  },
  {
   "cell_type": "raw",
   "id": "8b5135d6-869d-4717-8976-ff7f8daab103",
   "metadata": {
    "tags": []
   },
   "source": [
    "# line-dots colors                \n",
    "color_MOD = 'cyan'\n",
    "color_MYD = '#F2583E'\n",
    "color_MXD = 'magenta'\n",
    "color_THAAO = 'green'\n",
    "\n",
    "# general default settings\n",
    "x = ''\n",
    "y = ''\n",
    "xlabel = 'IWV MOD    [mm]' \n",
    "ylabel = 'IWV VESPA    [mm]'\n",
    "color = ''\n",
    "label = ''\n",
    "height = 500\n",
    "width = 600\n",
    "c = ''\n",
    "s = ''\n",
    "cmap = ''\n",
    "alpha = 10\n",
    "scale = 1\n",
    "size = 2\n",
    "legend = 'top_right'\n",
    "dynspread = False\n",
    "datashade = False  # datashading not available, può essere tolto\n",
    "grid=True\n",
    "title = ''\n",
    " \n",
    "\n",
    "splot_defaults = dict(x=x, y=x,\n",
    "                      xlabel=xlabel, ylabel=ylabel,\n",
    "                      color=color, label=label,\n",
    "                      height=height, width=width,\n",
    "                      s=s, c=c, \n",
    "                      cmap=cmap, \n",
    "                      alpha=alpha,\n",
    "                      scale=scale, size=size, \n",
    "                      legend=legend,\n",
    "                      dynspread=dynspread,datashade=datashade,\n",
    "                      grid=grid,\n",
    "                      title = title)\n",
    "\n",
    "# come cambiare i vaori delle keywords definite o aggiungerne altri\n",
    "# metodo 1\n",
    "print(splot_defaults)\n",
    "s_='diff_int'\n",
    "changed_plot_values = dict(s=s_)\n",
    "splot_defaults.update(changed_plot_values)\n",
    "print(f' metodo 1: {splot_defaults}')\n",
    "\n",
    "# metodo 2\n",
    "s_=''\n",
    "splot_defaults.update(s=s_)\n",
    "print(f' metodo 2: {splot_defaults}')\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6d630e92-2825-4add-9b11-e882431b4c42",
   "metadata": {
    "tags": []
   },
   "source": [
    "### **Basic scatter plots variables and visual kwargs options**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "24bef9a9-cd69-4ef7-98e7-b759a6adab12",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "splot_defaults.update(x='IWV_MODIS', y='IWV_THAAO', \n",
    "                      xlabel = 'IWV MOD    [mm]', \n",
    "                      ylabel = 'IWV VESPA    [mm]')\n",
    "splot_defaults.update(scale=10)\n",
    "print(f' splot_defaults: {splot_defaults}')\n",
    "\n",
    "splot_defaults.update(color=color_MOD, label='MOD')\n",
    "vespaVSmod = new_mod_vespa.hvplot.scatter(**splot_defaults)\n",
    "\n",
    "splot_defaults.update(color=color_MYD, label='MYD')\n",
    "vespaVSmyd = new_myd_vespa.hvplot.scatter(**splot_defaults)\n",
    "\n",
    "splot_defaults.update(color=color_MXD, label='MXD')\n",
    "vespaVSmxd = new_mxd_vespa.hvplot.scatter(**splot_defaults)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b37e4c95-e8eb-42d9-8dce-6f7a8fe902ec",
   "metadata": {
    "tags": []
   },
   "source": [
    "#### **Scatter plots: basic**"
   ]
  },
  {
   "cell_type": "raw",
   "id": "cea33a43-215b-46b2-a383-50ba249504ce",
   "metadata": {
    "tags": []
   },
   "source": [
    "# MOD e MYD separati     # ok: raw per test\n",
    "(vespaVSmod*slope*regress_mod_sk[0]*regress_mod_sk_text\n",
    "+ vespaVSmyd*slope*regress_myd_sk[0]*regress_myd_sk_text).cols(1)"
   ]
  },
  {
   "cell_type": "raw",
   "id": "fc811873-9493-4b5e-bd4f-3473c4e0705a",
   "metadata": {
    "tags": []
   },
   "source": [
    "# MXD: MOD e MYD insieme (stesso colore, stesso dataset mxd)    # ok: raw per test\n",
    "(vespaVSmxd*slope*regress_mxd_vespa_sk_text\n",
    "* regress_mod_sk[0]*regress_myd_sk[0]*regress_mxd_vespa_sk[0]\n",
    ").opts(xlabel='IWV MXD    [mm]', title=t_mxd_vespa, show_grid=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d526cfd3-c2b9-4169-a2db-773fd16f7817",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# MXD: MOD e MYD insieme (colori diversi, dataset mod e myd diversi)\n",
    "(vespaVSmod*vespaVSmyd*slope*regress_mxd_vespa_sk_text\n",
    "* regress_mod_sk[0]*regress_myd_sk[0]*regress_mxd_vespa_sk[0]\n",
    "  ).opts(xlabel='IWV MXD    [mm]', title=t_mxd_vespa, show_grid=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e8b3a600-1ff6-4062-b82d-3af352289a24",
   "metadata": {
    "tags": []
   },
   "source": [
    "### **Time differences (seconds) scatter plots variables and visual options**\n",
    "Scatter plots with the size of the dots that varies according to the <br>\n",
    "difference in seconds (diff_int) between the MODIS measurement and the <br>\n",
    "corresponding VESPA measurement selected within the +-deltatime range<br>\n",
    "\n",
    "**Da questi grafici si vede che la differenza di tempo in secondi tra le misure del  <br>\n",
    "punto MODIS e di VESPA non cambia nulla per i dati selezionati nell'intervallo  <br>\n",
    "deltatime=10m fissato.** <br>\n",
    "\n",
    "TODO: dovrebbe  essere anche il colore a variare in funzione della stessa <br>\n",
    "quantità, ma al momento sembra non funzionare, come anche la legenda<br>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7b344fda-c8b3-4a39-95f9-8bcdfe8fda29",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# Difference in time (seconds) between MODIS and VESPAdata points\n",
    "splot_defaults.update(c='diff_int', s='diff_int', \n",
    "                      scale=0.5, alpha=0.7)\n",
    "\n",
    "splot_defaults.update(color=color_MOD, label='MOD')\n",
    "vespaVSmod = new_mod_vespa.hvplot.scatter(**splot_defaults)\n",
    "\n",
    "splot_defaults.update(color=color_MYD, label='MYD')\n",
    "vespaVSmyd = new_myd_vespa.hvplot.scatter(**splot_defaults)\n",
    "\n",
    "splot_defaults.update(color=color_MXD, label='MXD')\n",
    "vespaVSmxd = new_mxd_vespa.hvplot.scatter(**splot_defaults)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3be2d1f6-a853-42ad-8c76-b18cb17a1b7c",
   "metadata": {
    "tags": []
   },
   "source": [
    "#### **Scatter plots: time difference (s)**\n",
    "Scatter plots with the size of the dots that varies according to the <br>\n",
    "difference in seconds (diff_int) between the MODIS measurement and the <br>\n",
    "corresponding VESPA measurement selected within the +-deltatime range<br>"
   ]
  },
  {
   "cell_type": "raw",
   "id": "8d3ce0c6-a57b-4bf0-ad4e-8981e151ca9a",
   "metadata": {
    "tags": []
   },
   "source": [
    "# MOD e MYD separati     # ok: raw per test\n",
    "(vespaVSmod*slope*regress_mod_sk[0]*regress_mod_sk_text\n",
    "+ vespaVSmyd*slope*regress_myd_sk[0]*regress_myd_sk_text).cols(1)"
   ]
  },
  {
   "cell_type": "raw",
   "id": "0a011c3a-bcb5-45b6-a0e7-9261a68645e4",
   "metadata": {
    "tags": []
   },
   "source": [
    "# MXD: MOD e MYD insieme (stesso colore, stesso dataset mxd)    # ok: raw per test\n",
    "(vespaVSmxd*slope*regress_mxd_vespa_sk_text\n",
    "*regress_mod_sk[0]*regress_myd_sk[0]*regress_mxd_vespa_sk[0])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d793cfd9-1713-4f01-b94f-182abd4a64f8",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# MXD: MOD e MYD insieme (colori diversi, dataset mod e myd diversi)\n",
    "(vespaVSmod*vespaVSmyd*slope*regress_mxd_vespa_sk_text\n",
    "*regress_mod_sk[0]*regress_myd_sk[0]*regress_mxd_vespa_sk[0])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "60e5c16d-85d7-49f8-88b3-bb90866034a0",
   "metadata": {},
   "source": [
    "##### **Same time differences (seconds) scatter plots with histograms**\n",
    "Scatter plots with the size of the dots that varies according to the <br>\n",
    "difference in seconds (diff_int) between the MODIS measurement and the <br>\n",
    "corresponding VESPA measurement selected within the +-deltatime range<br>\n",
    "with histograms <br>\n",
    "\n",
    "TODO: dovrebbe  essere anche il colore a variare in funzione della stessa <br>\n",
    "quantità, ma al momento sembra non funzionare, come anche la legenda<br>"
   ]
  },
  {
   "cell_type": "raw",
   "id": "e64e8553-af5e-48d6-8a7e-36cdbc9fcf66",
   "metadata": {
    "tags": []
   },
   "source": [
    "# Difference in time (seconds) between MODIS and VESPAdata points\n",
    "splot_defaults.update(c='diff_int', s='diff_int', \n",
    "                      scale=0.5, alpha=0.7)\n",
    "\n",
    "splot_defaults.update(color=color_MOD, label='MOD')\n",
    "vespaVSmod = new_mod_vespa.hvplot.scatter(**splot_defaults\n",
    "                                         ).hist(num_bins=21, dimension=['IWV_MODIS','IWV_THAAO'])\n",
    "\n",
    "splot_defaults.update(color=color_MYD, label='MYD')\n",
    "vespaVSmyd = new_myd_vespa.hvplot.scatter(**splot_defaults\n",
    "                                         ).hist(num_bins=21, dimension=['IWV_MODIS','IWV_THAAO'])\n",
    "\n",
    "splot_defaults.update(color=color_MXD, label='MXD')\n",
    "vespaVSmxd = new_mxd_vespa.hvplot.scatter(**splot_defaults\n",
    "                                         ).hist(num_bins=21, dimension=['IWV_MODIS','IWV_THAAO'])"
   ]
  },
  {
   "cell_type": "raw",
   "id": "e3e80e3e-5190-4870-a643-fe3d6b6c92b2",
   "metadata": {
    "tags": []
   },
   "source": [
    "# MOD e MYD separati     # ok: raw per test\n",
    "(vespaVSmod*slope*regress_mod_sk[0]*regress_mod_sk_text\n",
    "+ vespaVSmyd*slope*regress_myd_sk[0]*regress_myd_sk_text).cols(1)"
   ]
  },
  {
   "cell_type": "raw",
   "id": "81dded4c-c470-423a-96bc-93d3c3ef0209",
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    },
    "tags": []
   },
   "source": [
    "# MXD: MOD e MYD insieme (stesso colore, stesso dataset mxd)    # ok: raw per test\n",
    "(vespaVSmxd*slope*regress_mxd_vespa_sk_text\n",
    "*regress_mod_sk[0]*regress_myd_sk[0]*regress_mxd_vespa_sk[0])"
   ]
  },
  {
   "cell_type": "raw",
   "id": "38545290-ecfd-4f61-bfab-e5ea293841f9",
   "metadata": {
    "tags": []
   },
   "source": [
    "# MXD: MOD e MYD insieme (colori diversi, dataset mod e myd diversi)\n",
    "(vespaVSmod*vespaVSmyd*slope*regress_mxd_vespa_sk_text\n",
    "*regress_mod_sk[0]*regress_myd_sk[0]*regress_mxd_vespa_sk[0])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9277da7d",
   "metadata": {
    "tags": []
   },
   "source": [
    "##### Scatter dimensione punti su differenza secondi, OLD"
   ]
  },
  {
   "cell_type": "raw",
   "id": "9d4806af-38ca-41aa-842d-a9124ea2b1e0",
   "metadata": {
    "tags": []
   },
   "source": [
    "scale = 0.5\n",
    "alpha = 0.7\n",
    "vespaVSmod = new_mod_vespa.hvplot.scatter(x='IWV_MODIS', y='IWV_THAAO',\n",
    "                                      height=500, width=600, c='diff_int',\n",
    "                                      s='diff_int', scale=scale,\n",
    "                                      alpha=alpha,\n",
    "                                      xlabel='IWV MOD    [mm]', \n",
    "                                      ylabel='IWV THAAO    [mm]'\n",
    "                                     ).hist(num_bins=21,dimension=['IWV_MODIS','IWV_THAAO'])\n",
    "\n",
    "\n",
    "vespaVSmyd = new_myd_vespa.hvplot.scatter(x='IWV_MODIS', y='IWV_THAAO', \n",
    "                                      height=500, width=600, c='diff_int',\n",
    "                                      s='diff_int', scale=scale, \n",
    "                                      alpha=alpha,\n",
    "                                      xlabel='IWV MYD    [mm]', \n",
    "                                      ylabel='IWV THAAO    [mm]'\n",
    "                                     ).hist(num_bins=21,dimension=['IWV_MODIS','IWV_THAAO'])\n",
    "\n",
    "(vespaVSmod*slope*regress_mod_sk[0]*regress_mod_sk_text\n",
    "+ vespaVSmyd*slope*regress_myd_sk[0]*regress_myd_sk_text\n",
    "+ vespaVSmod*vespaVSmyd*slope*regress_mxd_vespa_sk_text\n",
    "* regress_mod_sk[0]*regress_myd_sk[0]*regress_mxd_vespa_sk[0]).cols(1)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4690d07f-583a-4626-af2a-60bf2fdcf708",
   "metadata": {},
   "source": [
    "### **Distance (km) scatter plots variables and visual options**\n",
    "Scatter plots with the size of the dots that varies according to the <br>\n",
    "distance in kilometers (diff_distance) between the MODIS point measurement and the <br>\n",
    "VESPA instrument location, # dots limited to those selected within the +-deltatime range<br>\n",
    "\n",
    "**Da questo si vede che la distanza tra il punto MODIS e la misura di VESPA nel punto fisso <br>\n",
    "del THAAO non cambia nulla per i punti nell'intervallo deltatime=10m fissato.** <br>\n",
    "\n",
    "TODO: dovrebbe  essere anche il colore a variare in funzione della stessa <br>\n",
    "quantità, ma al momento sembra non funzionare, come anche la legenda<br>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d03a0c91-14fd-421b-8f2e-9dc7e77509ec",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# Distance (km) between MODIS and VESPA data points\n",
    "splot_defaults.update(c='diff_distance', s='diff_distance', \n",
    "                      scale=5, alpha=0.7)\n",
    "\n",
    "splot_defaults.update(color=color_MOD, label='MOD')\n",
    "vespaVSmod = new_mod_vespa.hvplot.scatter(**splot_defaults)\n",
    "\n",
    "splot_defaults.update(color=color_MYD, label='MYD')\n",
    "vespaVSmyd = new_myd_vespa.hvplot.scatter(**splot_defaults)\n",
    "\n",
    "splot_defaults.update(color=color_MXD, label='MXD')\n",
    "vespaVSmxd = new_mxd_vespa.hvplot.scatter(**splot_defaults)"
   ]
  },
  {
   "cell_type": "raw",
   "id": "a2918485-aa80-459f-a4d4-228b7eb0afff",
   "metadata": {
    "jupyter": {
     "outputs_hidden": true
    },
    "tags": []
   },
   "source": [
    "# MOD e MYD separati     # ok: raw per test\n",
    "(vespaVSmod*slope*regress_mod_sk[0]*regress_mod_sk_text\n",
    "+ vespaVSmyd*slope*regress_myd_sk[0]*regress_myd_sk_text\n",
    " ).opts(width=400, height=300).cols(1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7e321d27-eb9d-4104-8fb1-a0aaf307ceba",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# MXD: MOD e MYD insieme (stesso colore, stesso dataset mxd)    # ok: raw per test\n",
    "# NON SO PERCHÉ: DOVREBBERO AVERE COLORI DIVERSI: DATASET DIVERSI MOD E MYD\n",
    "title = f'IWV MXD-VESPA              deltatime={deltatime}         # data points={num_mxd_data_p_Dt}'\n",
    "\n",
    "(vespaVSmxd*slope*regress_mxd_vespa_sk_text\n",
    " *regress_mod_sk[0]*regress_myd_sk[0]*regress_mxd_vespa_sk[0]\n",
    "#).opts(width=400, height=300\n",
    ").opts(title=title, show_grid=True, xlabel= 'IWV MXD    [mm]')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1d4f10e0-6cb8-4eef-9ba9-29d5cd260b9f",
   "metadata": {},
   "source": [
    "##### Scatter dimensione punti su differenza distanza geografica, OLD"
   ]
  },
  {
   "cell_type": "raw",
   "id": "66f32d89-0d19-4e91-8c57-5bb4d287182d",
   "metadata": {
    "tags": []
   },
   "source": [
    "# diff_distance non viene calcolata in closest_in_time() perché varia poco e variano poco i risultati, \n",
    "# quindi la cella va in errore. In caso di verifica reinserire la colonna 'diff_distance' nel dataframe  \n",
    "\n",
    "scale = 5\n",
    "alpha = 0.7\n",
    "vespaVSmod_scatter = new_mod_vespa.hvplot.scatter(x='IWV_MODIS', y='IWV_THAAO',\n",
    "                                              height=500, width=600, c='diff_distance',\n",
    "                                              s='diff_distance', scale=scale,\n",
    "                                              alpha=alpha,\n",
    "                                              xlabel='IWV MOD    [mm]', \n",
    "                                              ylabel='IWV THAAO    [mm]'\n",
    "                                             ).hist(num_bins=21,dimension=['IWV_MODIS','IWV_THAAO'])\n",
    "\n",
    "vespaVSmyd_scatter = new_myd_vespa.hvplot.scatter(x='IWV_MODIS', y='IWV_THAAO',\n",
    "                                              height=500, width=600, c='diff_distance',\n",
    "                                              s='diff_distance', scale=scale,\n",
    "                                              alpha=alpha,\n",
    "                                              xlabel='IWV MYD    [mm]',\n",
    "                                              ylabel='IWV THAAO    [mm]'\n",
    "                                             ).hist(num_bins=21,dimension=['IWV_MODIS','IWV_THAAO'])\n",
    "\n",
    "(vespaVSmod_scatter*slope*regress_mod_sk[0]*regress_mod_sk_text\n",
    "+ vespaVSmyd_scatter*slope*regress_myd_sk[0]*regress_myd_sk_text\n",
    "+ vespaVSmod_scatter*vespaVSmyd_scatter*slope*regress_mxd_vespa_sk_text\n",
    " *regress_mod_sk[0]*regress_myd_sk[0]*regress_mxd_vespa_sk[0]).opts(width=400, height=300).cols(1)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9b7ada25",
   "metadata": {},
   "source": [
    "### **Dipendenza da VZA**\n",
    "I grafici seguenti mostrano che non c'è dipendenza dal VZA per il valore di IWV <br>\n",
    "L'ipotesi più probabile è che il prodotto IWV MODIS di livello 2 sia già corretto <br> \n",
    "per l'air mass factor <br>"
   ]
  },
  {
   "cell_type": "raw",
   "id": "46ad1cfc-5f6c-4190-8829-3a109f9b4c51",
   "metadata": {
    "tags": []
   },
   "source": [
    "# diff_distance non viene calcolata in closest_in_time() perché varia poco e variano poco i risultati, \n",
    "# quindi la cella va in errore. In caso di verifica reinserire la colonna 'diff_distance' nel dataframe  \n",
    "scale = 5\n",
    "vespaVSmod_scatter = new_mod_vespa.hvplot.scatter(x='IWV_MODIS', y='IWV_THAAO',\n",
    "                                                  color=color_MOD,\n",
    "                                              height=height, width=width, \n",
    "                                              c='vza',\n",
    "                                              s='diff_distance', scale=scale,\n",
    "                                              xlabel='IWV MOD    [mm]',\n",
    "                                              ylabel='IWV THAAO    [mm]'\n",
    "                                             ).hist(num_bins=21,dimension=['IWV_MODIS','IWV_THAAO'])\n",
    "\n",
    "\n",
    "vespaVSmyd_scatter = new_myd_vespa.hvplot.scatter(x='IWV_MODIS', y='IWV_THAAO',\n",
    "                                                  color=color_MYD,\n",
    "                                              height=height, width=width, \n",
    "                                              c='vza',\n",
    "                                              s='diff_distance', scale=scale, \n",
    "                                              xlabel='IWV MYD    [mm]', \n",
    "                                              ylabel='IWV THAAO    [mm]'\n",
    "                                             ).hist(num_bins=21,dimension=['IWV_MODIS','IWV_THAAO'])\n",
    "\n",
    "(vespaVSmod_scatter*slope*regress_mod_sk[0]*regress_mod_sk_text\n",
    "+ vespaVSmyd_scatter*slope*regress_myd_sk[0]*regress_myd_sk_text\n",
    "+ vespaVSmod_scatter*vespaVSmyd_scatter*slope*regress_mxd_vespa_sk_text\n",
    " *regress_mod_sk[0]*regress_myd_sk[0]*regress_mxd_vespa_sk[0]).cols(1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "dbc4dd6f",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# ok s='diff_distance'\n",
    "height=400\n",
    "width=500\n",
    "scale = 3\n",
    "title = f'IWV MXD-VESPA              deltatime={deltatime}         # data points={num_mxd_data_p_Dt}'\n",
    "\n",
    "vespaVSmod_scatter_ok = new_mod_vespa.hvplot.scatter(x='IWV_MODIS', y='IWV_THAAO',\n",
    "                                              height=height, width=width, \n",
    "#                                              c='diff_distance',\n",
    "                                              c='vza',\n",
    "                                              s='diff_distance', \n",
    "                                              scale=scale,\n",
    "                                              xlabel='IWV MOD    [mm]',\n",
    "                                              ylabel='IWV VESPA    [mm]')\n",
    "\n",
    "vespaVSmyd_scatter_ok = new_myd_vespa.hvplot.scatter(x='IWV_MODIS', y='IWV_THAAO',\n",
    "                                              height=height, width=width, \n",
    "#                                              c='diff_distance',\n",
    "                                              c='vza',\n",
    "                                              s='diff_distance', \n",
    "                                              scale=scale,\n",
    "                                              xlabel='IWV MYD    [mm]', \n",
    "                                              ylabel='IWV VESPA    [mm]')\n",
    "\n",
    "(vespaVSmod_scatter_ok*vespaVSmyd_scatter_ok\n",
    "*slope\n",
    "*regress_mxd_vespa_sk_text\n",
    "#*regress_MOD_sk[0]*regress_MYD_sk[0]\n",
    "*regress_mxd_vespa_sk[0]).opts(title=title, show_grid=True, xlabel= 'IWV MXD    [mm]')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9cb57457-fd9b-4e94-ada2-2ff1c536adfb",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "scale = 5\n",
    "vespaVSmod_scatter = new_mod_vespa.hvplot.scatter(x='vza', y='diff_perc_IWV',\n",
    "                                              height=height, width=width, \n",
    "#                                              c='vza',\n",
    "#                                              s='diff_distance', scale=scale,\n",
    "                                              xlabel='vza',\n",
    "                                              ylabel='difference IWV THAAO'\n",
    "                                             ) #.hist(num_bins=21,dimension=['IWV_MODIS','IWV_THAAO'])\n",
    "\n",
    "\n",
    "vespaVSmyd_scatter = new_myd_vespa.hvplot.scatter(x='vza', y='diff_perc_IWV',\n",
    "                                              height=height, width=width, \n",
    "#                                              c='vza',\n",
    "#                                              s='diff_distance', scale=scale, \n",
    "                                              xlabel='vza', \n",
    "                                              ylabel='difference IWV THAAO'\n",
    "                                             ) #.hist(num_bins=21,dimension=['IWV_MODIS','IWV_THAAO'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "29094593-d215-4578-bd35-e004c803697a",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "#(vespaVSmod_scatter*slope*regress_mod_sk[0]*regress_mod_sk_text\n",
    "#+ vespaVSmyd_scatter*slope*regress_myd_sk[0]*regress_myd_sk_text\n",
    "#+ vespaVSmod_scatter*vespaVSmyd_scatter*slope*regress_mxd_vespa_sk_text\n",
    "# *regress_mod_sk[0]*regress_myd_sk[0]*regress_mxd_vespa_sk[0]).cols(1)\n",
    "(vespaVSmod_scatter*vespaVSmyd_scatter)#.cols(1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7a1d64e5",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# ok s='vza'\n",
    "height=400\n",
    "width=500\n",
    "scale = 1# 0.5\n",
    "title_mxd = f'IWV MXD-VESPA              deltatime={deltatime}         # data points={num_mxd_data_p_Dt}'\n",
    "\n",
    "vespaVSmod_scatter_ok = new_mod_vespa.hvplot.scatter(x='IWV_MODIS', y='IWV_THAAO',\n",
    "                                              height=height, width=width, \n",
    "                                              c='vza',\n",
    "                                              s='vza', scale=scale,\n",
    "                                              xlabel='IWV MOD    [mm]',\n",
    "                                              ylabel='IWV VESPA    [mm]')\n",
    "\n",
    "vespaVSmyd_scatter_ok = new_myd_vespa.hvplot.scatter(x='IWV_MODIS', y='IWV_THAAO',\n",
    "                                              height=height, width=width, \n",
    "                                              c='vza',\n",
    "                                              s='vza', scale=scale, \n",
    "                                              xlabel='IWV MYD    [mm]', \n",
    "                                              ylabel='IWV VESPA    [mm]')\n",
    "\n",
    "(vespaVSmod_scatter_ok*vespaVSmyd_scatter_ok\n",
    "*slope\n",
    "*regress_mxd_vespa_sk_text\n",
    "#*regress_MOD_sk[0]*regress_MYD_sk[0]\n",
    "*regress_mxd_vespa_sk[0]).opts(title=title_mxd, show_grid=True, xlabel= 'IWV MXD    [mm]')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "883ab440",
   "metadata": {},
   "source": [
    "## **Stagionalità**"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "79e770f2-03c2-4329-8bbc-303b41846fa9",
   "metadata": {
    "tags": []
   },
   "source": [
    "### **Scatter plot GroupBy month 1 (con slider)**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "447531c3-d3da-4541-b8f7-ce80fb697b6b",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "height=300\n",
    "width=600\n",
    "xlim = [0,2]\n",
    "ylim = [0,2]\n",
    "xlim = [0,25]\n",
    "ylim = [0,25]\n",
    "#x_offset=xlim[1]-30\n",
    "#y_offset=ylim[0]+5\n",
    "scale = 5\n",
    "\n",
    "regress_mod_vespa_s_text = regression_eq_text(new_mod_vespa.IWV_MODIS,\n",
    "                                        new_mod_vespa.IWV_THAAO, regress_mod_sk,\n",
    "                                        x_offset=xlim[0]-0.5,y_offset=ylim[0]+0.5)\n",
    "regress_myd_vespa_s_text = regression_eq_text(new_myd_vespa.IWV_MODIS,\n",
    "                                        new_myd_vespa.IWV_THAAO, regress_myd_sk,\n",
    "                                        x_offset=xlim[0]-0.5,y_offset=ylim[0]+0.5)\n",
    "regress_mxd_vespa_s_text = regression_eq_text(new_mxd_vespa.IWV_MODIS,\n",
    "                                        new_mxd_vespa.IWV_THAAO, regress_mxd_vespa_sk,\n",
    "                                        x_offset=xlim[0]-0.5,y_offset=ylim[0]+0.5)\n",
    "\n",
    "vespaVSmod_scatter = new_mod_vespa.hvplot.scatter(x='IWV_MODIS', y='IWV_THAAO', \n",
    "                                              height=height, width=width, \n",
    "                                              c='diff_distance', \n",
    "                                              xlim=xlim, ylim=ylim, \n",
    "#                                              s='diff_distance', \n",
    "                                              scale=scale, \n",
    "                                              xlabel='IWV MOD    [mm]', \n",
    "                                              ylabel='IWV THAAO    [mm]',\n",
    "                                              groupby=['Time.month'])\n",
    "\n",
    "vespaVSmyd_scatter = new_myd_vespa.hvplot.scatter(x='IWV_MODIS', y='IWV_THAAO', \n",
    "                                              height=height, width=width, \n",
    "                                              c='diff_distance',\n",
    "                                              xlim=xlim, ylim=ylim, \n",
    "#                                              s='diff_distance', \n",
    "                                              scale=scale, \n",
    "                                              xlabel='IWV MYD    [mm]',\n",
    "                                              ylabel='IWV THAAO    [mm]',\n",
    "                                              groupby=['Time.month'])"
   ]
  },
  {
   "cell_type": "raw",
   "id": "560e792a-e4b1-41c2-800b-930601062bac",
   "metadata": {
    "tags": []
   },
   "source": [
    "# MOD e MYD separati     # ok: raw per test\n",
    "(vespaVSmod_scatter*slope*regress_mod_sk[0]*regress_mod_vespa_s_text\n",
    "+ vespaVSmyd_scatter*slope*regress_myd_sk[0]*regress_myd_vespa_s_text).cols(1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "cffa1650-be12-4fb6-a53a-d7b959730bb6",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# MOD e MYD insieme \n",
    "(vespaVSmod_scatter*vespaVSmyd_scatter*slope*regress_mxd_vespa_s_text\n",
    "*regress_mod_sk[0]*regress_myd_sk[0]*regress_mxd_vespa_sk[0]\n",
    ").opts(title=title_mxd, xlabel= 'IWV MXD    [mm]', show_grid=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "72f4706d",
   "metadata": {},
   "source": [
    "### GroupBy month"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0ed407e8-8687-418f-90be-e31e9a629324",
   "metadata": {
    "tags": []
   },
   "source": [
    "### **Scatter plot GroupBy month 2 (con slider)**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1a532c07-209e-4597-bc8d-adb8d866c035",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "height=300\n",
    "width=600\n",
    "\n",
    "vespaVSmod_scatter = new_mod_vespa.hvplot.scatter(x='IWV_MODIS', y='IWV_THAAO',\n",
    "                                              height=height, width=width, \n",
    "                                              xlim=xlim, ylim=ylim,\n",
    "                                              scale=0.3, xlabel='IWV MOD    [mm]',\n",
    "                                              ylabel='IWV THAAO    [mm]',\n",
    "                                              by=['Time.month'], groupby=['Time.month'])\n",
    "\n",
    "vespaVSmyd_scatter = new_myd_vespa.hvplot.scatter(x='IWV_MODIS', y='IWV_THAAO',\n",
    "                                              height=height, width=width, \n",
    "                                              xlim=xlim, ylim=ylim,\n",
    "                                              scale=0.3, xlabel='IWV MYD    [mm]',\n",
    "                                              ylabel='IWV THAAO    [mm]',\n",
    "                                              by=['Time.month'], groupby=['Time.month'])"
   ]
  },
  {
   "cell_type": "raw",
   "id": "5251009b-a7ef-4003-adc0-62249caada54",
   "metadata": {
    "tags": []
   },
   "source": [
    "# MOD e MYD separati     # ok: raw per test\n",
    "(vespaVSmod_scatter*slope*regress_mod_sk[0]*regress_mod_vespa_s_text\n",
    "+ vespaVSmyd_scatter*slope*regress_myd_sk[0]*regress_myd_vespa_s_text).cols(1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7934cc2b-c020-4dc9-9e52-cbe9856dd4e0",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "(vespaVSmod_scatter*vespaVSmyd_scatter*slope*regress_mod_sk[0]*regress_myd_sk[0]\n",
    "*regress_mxd_vespa_sk[0]*regress_mxd_vespa_s_text\n",
    ").opts(title='MXD vs VESPA IWV seasonality by month', \n",
    "#         xlabel = 'IWV MXD (mm)', ylabel='IWV THAAO (mm)')).cols(1)\n",
    "       xlabel = 'IWV MXD    [mm]', ylabel='IWV THAAO    [mm]', show_grid=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "908e59a2-d3fc-4bc9-be29-ef9012d452b4",
   "metadata": {},
   "source": [
    "### **Scatter plot GroupBy month (animation)**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "45c68b54-81c6-4060-a079-9f44bdb32068",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "vespaVSmod_scatter = new_mod_vespa.hvplot.scatter(x='IWV_MODIS', y='IWV_THAAO',\n",
    "                                              height=300,\n",
    "                                              xlim=xlim, ylim=ylim,\n",
    "                                              scale=0.3, \n",
    "                                              xlabel='IWV MOD    [mm]', ylabel='IWV THAAO    [mm]',\n",
    "                                              by=['Time.month'], groupby=['Time.month'])\n",
    "\n",
    "vespaVSmyd_scatter = new_myd_vespa.hvplot.scatter(x='IWV_MODIS', y='IWV_THAAO',\n",
    "                                              height=300, \n",
    "                                              xlim=xlim, ylim=ylim,\n",
    "                                              scale=0.3,\n",
    "                                              xlabel='IWV MYD    [mm]', ylabel='IWV THAAO    [mm]',\n",
    "                                              by=['Time.month'], groupby=['Time.month'])"
   ]
  },
  {
   "cell_type": "raw",
   "id": "47d870e6-6f7e-4315-b0fc-c6354b602cc8",
   "metadata": {
    "tags": []
   },
   "source": [
    "# MODe MYD separati : ok per test \n",
    "pn.panel((vespaVSmod_scatter*slope*regress_mod_sk[0]*regress_mod_vespa_s_text\n",
    "+ vespaVSmyd_scatter*slope*regress_myd_sk[0]*regress_myd_vespa_s_text\n",
    ").cols(1),\n",
    "#         xlabel = 'IWV MXD (mm)', ylabel='IWV THAAO (mm)')).cols(1),\n",
    "       widget_type='scrubber', widget_location='top')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1a9d1d23-35cc-4283-ab0d-0dfd6d25b81e",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "pn.panel((vespaVSmod_scatter*vespaVSmyd_scatter*slope*regress_mod_sk[0]*regress_myd_sk[0]\n",
    "*regress_mxd_vespa_sk[0]*regress_mxd_vespa_s_text\n",
    ").opts(title='MXD vs VESPA IWV seasonality animation', \n",
    "       xlabel = 'IWV MXD    [mm]', ylabel='IWV THAAO    [mm]', show_grid=True),\n",
    "       widget_type='scrubber', widget_location='top')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ea89e210",
   "metadata": {},
   "source": [
    "### **Scatter plot By month**"
   ]
  },
  {
   "cell_type": "raw",
   "id": "3f6500f5-a711-454b-ab04-c82048faf4ab",
   "metadata": {
    "tags": []
   },
   "source": [
    "# MOD e MYD separati     # ok: raw per test\n",
    "#cmap = 'cet_CET_C5'\n",
    "height=349\n",
    "width=500\n",
    "vespaVSmod_scatter = new_mod_vespa.hvplot.scatter(x='IWV_MODIS', y='IWV_THAAO',\n",
    "                                               height=height, width=width,\n",
    "                                              cmap=cmap,\n",
    "                                              scale=0.3, \n",
    "                                              xlabel='IWV MOD    [mm]', ylabel='IWV THAAO    [mm]',\n",
    "                                              by=['Time.month'])\n",
    "\n",
    "vespaVSmyd_scatter = new_myd_vespa.hvplot.scatter(x='IWV_MODIS', y='IWV_THAAO',\n",
    "                                              height=height, width=width,\n",
    "                                              cmap=cmap,\n",
    "                                              scale=0.3, \n",
    "                                              xlabel='IWV MOD    [mm]', ylabel='IWV THAAO    [mm]',\n",
    "                                              by=['Time.month'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3ae98f20-5a4e-421c-af4d-a0fc14b530a6",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "(vespaVSmod_scatter*slope*regress_mod_sk[0]*regress_mod_vespa_s_text\n",
    "+ vespaVSmyd_scatter*slope*regress_myd_sk[0]*regress_myd_vespa_s_text\n",
    ").cols(1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b627525f",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "# ok stagionalità\n",
    "#cmap = 'hsv'\n",
    "height=370\n",
    "width=500\n",
    "vespaVSmod_scatter = new_mod_vespa.hvplot.scatter(x='IWV_MODIS', y='IWV_THAAO',\n",
    "                                               height=height, width=width,\n",
    "#                                              cmap=cmap,\n",
    "                                              scale=0.3, \n",
    "                                              xlabel='IWV MOD    [mm]', ylabel='IWV VESPA    [mm]',\n",
    "                                              by=['Time.month'])\n",
    "\n",
    "vespaVSmyd_scatter = new_myd_vespa.hvplot.scatter(x='IWV_MODIS', y='IWV_THAAO',\n",
    "                                              height=height, width=width,\n",
    "#                                              cmap=cmap,\n",
    "                                              scale=0.3, \n",
    "                                              xlabel='IWV MOD    [mm]', ylabel='IWV VESPA    [mm]',\n",
    "                                              by=['Time.month'])\n",
    "\n",
    "\n",
    "(vespaVSmod_scatter*vespaVSmyd_scatter\n",
    "*slope\n",
    "#*regress_MOD_sk[0]*regress_MYD_sk[0]\n",
    "*regress_mxd_vespa_sk[0]*regress_mxd_vespa_sk_text\n",
    ").opts(title='IWV seasonality', \n",
    "       xlabel = 'IWV MXD    [mm]', ylabel='IWV THAAO    [mm]', show_grid=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "990e7f97",
   "metadata": {},
   "source": [
    "## **Diff IWV MXD-THAAO vs IWV THAAO** "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d828e93a",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "height=370\n",
    "width=800"
   ]
  },
  {
   "cell_type": "raw",
   "id": "376a0ea4-8c20-4067-b6f0-d83cb97edb02",
   "metadata": {},
   "source": [
    "new_mod_vespa.head(1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "038e431f-4755-45dd-8c71-dab6c1b40fbd",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "Diff_iwv_vespaVSmod_scatter = new_mod_vespa.hvplot.scatter(x='IWV_THAAO',\n",
    "                                                   y='diff_IWV', by=['Time.month'],\n",
    "                                                   xlabel='IWV VESPA    [mm]',\n",
    "                                                   ylabel='Diff MOD-VESPA    [mm]',\n",
    "                                                   height=height, width=width, \n",
    "                                                   scale=scale,grid=True)\n",
    "\n",
    "Diff_iwv_vespaVSmyd_scatter = new_myd_vespa.hvplot.scatter(x='IWV_THAAO',\n",
    "                                                   y='diff_IWV', by=['Time.month'],\n",
    "                                                   xlabel='IWV VESPA    [mm]',\n",
    "                                                   ylabel='Diff MYD-VESPA    [mm]',\n",
    "                                                   height=height, width=width, \n",
    "                                                   scale=scale,grid=True)\n",
    "\n",
    "(Diff_iwv_vespaVSmod_scatter\n",
    "+ Diff_iwv_vespaVSmyd_scatter\n",
    "+ (Diff_iwv_vespaVSmod_scatter*Diff_iwv_vespaVSmyd_scatter\n",
    "  ).opts(ylabel='Diff MXD-VESPA    [mm]')).cols(1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d823a390",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "Diff_iwv_vespaVSmod_scatter = new_mod_vespa.hvplot.scatter(x='IWV_THAAO',\n",
    "#                                                   y='diff_IWV', by=['time.month'],\n",
    "#                                                   xlabel='IWV VESPA    [cm]',\n",
    "#                                                   ylabel='Diff MOD-VESPA    [cm]',\n",
    "                                                   y='diff_IWV', by=['Time.month'],\n",
    "                                                   xlabel='IWV VESPA    [mm]',\n",
    "                                                   ylabel='Diff MOD-VESPA    [mm]',\n",
    "                                                   height=height, width=width, \n",
    "                                                   scale=scale,grid=True)\n",
    "\n",
    "Diff_iwv_vespaVSmyd_scatter = new_myd_vespa.hvplot.scatter(x='IWV_THAAO',\n",
    "#                                                   y='diff_IWV', by=['time.month'],\n",
    "#                                                   xlabel='IWV VESPA    [cm]',\n",
    "#                                                   ylabel='Diff MYD-VESPA    [cm]',\n",
    "                                                   y='diff_IWV', by=['Time.month'],\n",
    "                                                   xlabel='IWV VESPA    [mm]',\n",
    "                                                   ylabel='Diff MYD-VESPA    [mm]',\n",
    "                                                   height=height, width=width, \n",
    "                                                   scale=scale,grid=True)\n",
    "\n",
    "(Diff_iwv_vespaVSmod_scatter\n",
    "+ Diff_iwv_vespaVSmyd_scatter\n",
    "+ (Diff_iwv_vespaVSmod_scatter*Diff_iwv_vespaVSmyd_scatter\n",
    "  ).opts(ylabel='Diff MXD-VESPA    [mm]')).cols(1)"
   ]
  },
  {
   "cell_type": "raw",
   "id": "9038f38f",
   "metadata": {},
   "source": [
    "pippo = (Diff_iwv_Vs_MOD_scatter*Diff_iwv_vespaVSmyd_scatter).opts(ylabel='Diff MXD-THAAO (cm)')\n",
    "pippo = pn.panel(pippo)\n",
    "pippo.save('pippo.html')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e062358d",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "title = f'IWV MXD-VESPA              deltatime={deltatime}             # data points={num_mxd_data_p_Dt}'\n",
    "((Diff_iwv_vespaVSmod_scatter*Diff_iwv_vespaVSmyd_scatter).\n",
    " opts(title=title, ylabel='Diff MXD-VESPA    [mm]'))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d70dfbf2",
   "metadata": {},
   "source": [
    "### Diff % IWV MXD-THAAO vs IWV THAAO "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1d38b28e-7ef6-4686-bf9b-d09bfdaa517b",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "Diffperc_iwv_vespaVSmod_scatter = new_mod_vespa.hvplot.scatter(x='IWV_THAAO',\n",
    "                                                   y='diff_perc_IWV', by=['Time.month'],\n",
    "                                                   xlabel='IWV VESPA    [mm]',\n",
    "                                                   ylabel='100*(MOD-VESPA)/VESPA',\n",
    "                                                   height=height, width=width, \n",
    "                                                   scale=scale,grid=True)\n",
    "\n",
    "Diffperc_iwv_vespaVSmyd_scatter = new_myd_vespa.hvplot.scatter(x='IWV_THAAO',\n",
    "                                                   y='diff_perc_IWV', by=['Time.month'],\n",
    "                                                   xlabel='IWV VESPA    [mm]',\n",
    "                                                   ylabel='100*(MYD-VESPA)/VESPA',\n",
    "                                                   height=height, width=width, \n",
    "                                                   scale=scale,grid=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "433fdfa5-44f1-4b63-a83f-b01e0cc20faa",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "(Diffperc_iwv_vespaVSmod_scatter\n",
    "+ Diffperc_iwv_vespaVSmyd_scatter).cols(1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6444b4f3-71a7-483e-9778-a52db54ff396",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "(Diffperc_iwv_vespaVSmod_scatter*Diffperc_iwv_vespaVSmyd_scatter\n",
    ").opts(ylabel='100*(MXD-VESPA)/VESPA')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "8c1ae4e0",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "title = f'IWV 100*(MXD-VESPA)/VESPA              deltatime={deltatime}             # data points={num_mxd_data_p_Dt}'      # data points={number_MXD_data_points_deltatime}'\n",
    "((Diffperc_iwv_vespaVSmod_scatter*Diffperc_iwv_vespaVSmyd_scatter)\n",
    " .opts(title=title, ylabel='100*(MXD-VESPA)/VESPA'))"
   ]
  },
  {
   "cell_type": "raw",
   "id": "f26ce18d",
   "metadata": {},
   "source": [
    "type(read_hatpro_year('/home/luca/WORK/DATA/Groundbased/THULE/HATPRO/IWV_2017_HATPRO.txt'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "71502054-e9c8-4656-977a-506c8e87f542",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "",
   "language": "python",
   "name": ""
  },
  "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.16"
  },
  "toc-autonumbering": true,
  "toc-showmarkdowntxt": false,
  "toc-showtags": true
 },
 "nbformat": 4,
 "nbformat_minor": 5
}