Merge pull request #3 from LorenzoStucchi/full-auto

LorenzoStucchi · web-flow · commit bf65b85d2a2d · 2021-11-04T18:28:18.000+01:00
Full auto methodology
diff --git a/README.md b/README.md
@@ -1,9 +1,9 @@
 # Automatization of the SEN-et SNAP plugin for multi-temporal analysis
 ## Description
 
-The [SEN-ET project](https://www.esa-sen4et.org/) has been founded for ESA with the main object to "develop an optimal methodology for estimating evapotranspiration at fine (tens of meters) spatial scale, based on synergistic use of Sentinel 2 and Sentinel 3 satellites’ observations" as reported by the website of the project.
+The [SEN-ET project](https://www.esa-sen4et.org/) has been founded for ESA with the main object to *"develop an optimal methodology for estimating evapotranspiration at fine (tens of meters) spatial scale, based on synergistic use of Sentinel 2 and Sentinel 3 satellites’ observations"* as reported by the website of the project.
 
-The plugin is mainly thought for use inside SNAP, so a manual operation where the user has to make the different operations manually. However, also the python script for the operation have been released.
+The plugin is mainly thought for use inside SNAP, so a manual operation where the user has to make the different operations manually. However, also the python scripts for the operation have been released.
 The scope of this repository is to create the possibility to automatize the process, allowing the computation of multiple days without the manual input of the operator. So, it is possible to compute monthly means on an area of interest.
 
 The original schema for the plugin is represented in the image below:
@@ -12,14 +12,36 @@ The original schema for the plugin is represented in the image below:
 
 The automation process reduces it to one step and presents this schema:
 
-![New short sen-et](assets/new_one_step.png)
+![New short sen-et](assets/new_full.png)
 
 ## Installation
 It is needed to have installed SNAP (with Python 3.6) and the [SEN-ET plugin](https://www.esa-sen4et.org/static/media/Sen-ET-plugin-v1.0.1.b41ae6c8.zip) as indicated in the [official plugin documentation](https://www.esa-sen4et.org/static/media/sen-et-user-manual-v1.1.0.5d1ac526.pdf) in Section 3.2. In case of troubleshooting, see Section 3.4.
 
-The test that everything is working will be the download of the ERA5 Data in the next step.
+The suggested version of python to be used for Windows is the 3.6.8, it should be download before the installation of SNAP. The installer could installed [here](https://www.python.org/ftp/python/3.6.8/python-3.6.8-amd64.exe).
 
-It is also needed for Windows computers to install the Git Bash terminal provided with the installation of [Git](https://git-scm.com/downloads).
+At the first start, it is needed to run an operator to complete the bundle installation as indicated in the documentation:
+
+*The first time any of the Sen-ET operators is run the user will receive a SNAP - Warning message about path not existing. After clicking OK another dialog will appear asking user if they want to proceed with bundle download/installation (Figure 3.4). After clicking "Yes" the download and installation of the plugin bundle (the self-contained Python environment) will happen automatically. When the bundle installation finishes, a message will notify of the successful installation of the bundle. The Sen-ET SNAP plugin is installed and ready for use.*
+
+The installation and the test that everything is working will be done in the step of the download of the ERA5 Data. This will be later explained.
+
+For Windows computers, it is also needed to install the Git Bash terminal provided with the installation of [Git](https://git-scm.com/downloads).
+
+Download the code in the terminal with
+```
+git clone https://github.com/LorenzoStucchi/SEN-ET_automatization.git
+```
+
+### Python for generating scripts
+For the automatization of the project, it is needed to have python installed on the computer. Also, these libraries are needed:
+  - xarray
+  - netcdf4
+  - h5netcdf
+
+It is possible and recommended to use conda to create the environment with all the dependencies with the command:
+```
+conda env create -f env.yml
+```
 
 ## Data download
 ### ECMWF Reanalysis v5 (ERA5) Data
@@ -43,7 +65,7 @@ The images request needs to be of the product type:
 That option could be added in the filters, zoom to the area of interest and drawn a rectangular in the area. Open the filter panel and add in the sensing date, the timespan of interest for your analysis, and the desired timespan. Using the preview, it is also possible to check the cloud coverage visually.
 
 ### Folder structure
-Once the data are downloaded, it is needed to organize them in a specific way. An analysis of the data for one month will produce products for around 500 GB. The main folder that could be named as prefered, in the schema below `MAIN_FOLDER`. One folder for each type of data called `era`, `S2` and `S3`. Inside the `S2` and `S3` folders, a folder for each of the days of the images, the folder should be named in this format `YYYY_MM_DD`.
+Once the data are downloaded, it is needed to organize them in a specific way. An analysis of the data for one month will produce products for around 500 GB. The main folder that could be named as preferred, in the schema below `MAIN_FOLDER`. One folder for each type of data is called `era`, `S2` and `S3`. Inside the `S2` and `S3` folders, a folder for each of the days of the images, the folder should be named in this format `YYYY_MM_DD`.
 
 Once the structure has been created, the `era5.nc` and the images should be moved to the correct folder. Also, it is needed to unzip the images.
 As example:
@@ -88,25 +110,31 @@ MAIN_FOLDER
 Set the parameters in the file [`parameters.json`](input/parameters.json):
 - Set the coordinates of the Area Of Interest (AOI), which could be the portion of the Sentinel 2 images or a bigger area; it is used to cut the Sentinel 3 images and improve the speed of the computation.
 - The absolute path indicated by `MAIN_FOLDER` in the previous schema should be inserted in the `general_path` variable. In the example, the `MAIN_FOLDER` absolute path is `"D:\\TEMP\\"`.
-- The path into the `setet_folder` variable. The `setet_folder` is something like `"C:\\Users\\user\\.snap\\auxdata\\sen-et-conda-Win64"`.
+- The path into the `senet_folder` variable. The `senet_folder` is something like `"C:\\Users\\user\\.snap\\auxdata\\sen-et-conda-Win64"`.
 - It is also possible to modify all the default computational parameters, as explained in Section 3.3 of the official guide. The *only* parameter that should be modified is the [timezone of the area](input/parameters.json#L12).
 
-Finally, copy the absolute path of all the Sentinel Images in the file [`s2_paths.txt`](input/s2_paths.txt) and [`s3_paths.txt`](input/s3_paths.txt). The path should point to the file `MTD_MSIL2A.xml` for S2 and `xfdumanifest.xml` for S3.
+Finally, copy the absolute path of all the Sentinel Images in the files [`s2_paths.txt`](input/s2_paths.txt) and [`s3_paths.txt`](input/s3_paths.txt). The path should point to the folders `.SAFE` for S2 and `.SEN3` for S3.
 
 ## Run the code
-In a Bash terminal run: 
+### Windows
+In conda terminal run: 
 ```
-sh main.sh
+conda activate senet-auto
+main.sh
 ```
-After the computation of the Sentinel 3 images in the graph processing tool, it is necessary to insert the time of acquisition for the Sentinel 3 images. So, it is needed to open in SNAP the Sentinel image cut on the area of interest and report the time in the file `days_ADDTIME_s3.txt` and rename it as `days_time_s3.txt`.
-The structure of the file should be: `YYYY_MM_DD HH:MM`
+This will open a bash Terminal with all the dependencies correctly installed. 
 
-The script produces a message when the Sentinel 3 images are available for the computation and requests to confirm that the file `days_time_s3.txt` has been created.
+### Unix system (linux/MacOS)
+In a terminal run: 
+```
+conda activate senet-auto
+sh main.sh
+```
 
 ## Correction of errors
 In the [official code repository](https://github.com/DHI-GRAS/sen-et-snap-scripts) the graph presents some errors, due to some changes in the new versions of SNAP. To fix the errors the new version of the graph is present in the folder [graph](graph/).
 
 ## Authors
 The official code of SNAP is released with [GNU General Public License v3.0](https://github.com/DHI-GRAS/sen-et-snap-scripts/blob/master/LICENSE) by the original authors and as reported in the [plugin code repository](https://github.com/DHI-GRAS/senEtSnapSta).
 
-The code for the automatization of the process is also realeased by [Lorenzo Stucchi](https://github.com/LorenzoStucchi) with [GNU General Public License v3.0](LICENSE).
+The code for the automatization of the process is also released by [Lorenzo Stucchi](https://github.com/LorenzoStucchi) with [GNU General Public License v3.0](LICENSE).
diff --git a/assets/new.png b/assets/new.png
diff --git a/assets/new_full.png b/assets/new_full.png
diff --git a/assets/new_one_step.png b/assets/new_one_step.png
diff --git a/assets/old.png b/assets/old.png
diff --git a/env.yml b/env.yml
@@ -0,0 +1,9 @@
+name: senet-auto
+channels:
+  - conda-forge
+  - defaults
+dependencies:
+  - python=3.9
+  - xarray
+  - netcdf4
+  - h5netcdf
diff --git a/graph/sentinel_3_cut_time.xml b/graph/sentinel_3_cut_time.xml
@@ -24,6 +24,33 @@
       <copyMetadata>true</copyMetadata>
     </parameters>
   </node>
+  <node id="Subset_remove_bands">
+    <operator>Subset</operator>
+    <sources>
+      <sourceProduct refid="Subset"/>
+    </sources>
+    <parameters class="com.bc.ceres.binding.dom.XppDomElement">
+      <sourceBands>LST</sourceBands>
+      <region/>
+      <referenceBand/>
+      <geoRegion/>
+      <subSamplingX>1</subSamplingX>
+      <subSamplingY>1</subSamplingY>
+      <fullSwath>false</fullSwath>
+      <tiePointGridNames/>
+      <copyMetadata>true</copyMetadata>
+    </parameters>
+  </node>
+  <node id="Write_netcdf">
+    <operator>Write</operator>
+    <sources>
+      <sourceProduct refid="Subset_remove_bands"/>
+    </sources>
+    <parameters class="com.bc.ceres.binding.dom.XppDomElement">
+      <file>!OUTPUT_Sentinel-3_time!</file>
+      <formatName>NetCDF4-BEAM</formatName>
+    </parameters>
+  </node>
   <node id="Write">
     <operator>Write</operator>
     <sources>
@@ -42,6 +69,12 @@
     <node id="Subset">
       <displayPosition x="234.0" y="130.0"/>
     </node>
+    <node id="Subset_remove_bands">
+      <displayPosition x="233.0" y="224.0"/>
+    </node>
+    <node id="Write_netcdf">
+      <displayPosition x="458.0" y="226.0"/>
+    </node>
     <node id="Write">
             <displayPosition x="455.0" y="130.0"/>
     </node>
diff --git a/input/s2_paths.txt b/input/s2_paths.txt
@@ -1 +1 @@
-D:\TEMP\S2\2019_08_06\S2A_MSIL2A_20190806T095031_N0213_R079_T33TWF_20190806T114240.SAFE\MTD_MSIL2A.xml
+D:\TEMP\S2\2019_08_06\S2A_MSIL2A_20190806T095031_N0213_R079_T33TWF_20190806T114240.SAFE
diff --git a/input/s3_paths.txt b/input/s3_paths.txt
@@ -1 +1 @@
-D:\TEMP\S3\2019_08_03\S3A_SL_2_LST____20190803T083233_20190803T101332_20190804T153249_6059_047_335______LN2_O_NT_003.SEN3\xfdumanifest.xml
+D:\TEMP\S3\2019_08_03\S3A_SL_2_LST____20190803T083233_20190803T101332_20190804T153249_6059_047_335______LN2_O_NT_003.SEN3
diff --git a/main.sh b/main.sh
@@ -1,30 +1,10 @@
 #!/bin/bash
 
-python py_script/1_create_script_graph.py
+python py_scripts/1_create_script_graph.py
 sh output/1_script_gpt_s3.sh
-echo "      ===================================================================="
-echo "      ===================================================================="
-echo "      Open SNAP and set the hour of s3 images in the file days_time_s3.txt"
-echo "      ===================================================================="
-echo "      ===================================================================="
-sleep 30
+python py_scripts/2_read_times_s3.py
 sh output/2_script_gpt_s2.sh
-read -r -p "Do you set the S3 time in the file days_time_s3.txt? [Y/n] " input
-if [[ $input == Y ]] || [[ $input == y ]]
-then
-    python py_script/2_create_script_python.py
-    sh output/3_script_S2.sh
-    sh output/4_script_S3.sh
-    echo "Computation ended sucessfully"
-else
-    read -r -p "Please, set the hour and minute of the Sentinel 3 images in the file days_time_s3.txt and then press Y: " input
-    if [[ $input == Y ]] || [[ $input == y ]]
-    then
-        python py_script/2_create_script_python.py
-        sh output/3_script_S2.sh
-        sh output/4_script_S3.sh
-        echo "Computation ended sucessfully"
-    else
-        echo "Wrong value chose, the script ended without computation"
-    fi
-fi
+python py_scripts/3_create_script_python.py
+sh output/3_script_S2.sh
+sh output/4_script_S3.sh
+read -r -p "Computation ended sucessfully! Press any key to close the panel" input
diff --git a/py_scripts/1_create_script_graph.py b/py_scripts/1_create_script_graph.py
@@ -10,11 +10,11 @@
 path_zen = "graph\\sun_zenith_sentinel_2_pre_processing.xml"
 path_LC = "graph\\add_landcover.xml"
 path_ele = "graph\\add_elevation.xml"
-path_cut = "graph\\sentinel_3_cut.xml"
+path_cut = "graph\\sentinel_3_cut_time.xml"
 path_s3_pro = "graph\\sentinel_3_pre_processing_fix.xml"
 path_graph = "output\\graph\\"
 path_date_s2 = "output\\days_s2.txt"
-path_date_s3 = "output\\days_ADDTIME_s3.txt"
+path_date_s3 = "output\\days_s3.txt"
 path_script_s3 = "output\\1_script_gpt_s3.sh"
 path_script_s2 = "output\\2_script_gpt_s2.sh"
 
@@ -67,15 +67,18 @@
 date_s3 = ""
 for line in file_path_s3:
     s3_image_path = line.strip()
-    path = "\\".join(s3_image_path.split("\\")[:-2])
-    date = " ".join(s3_image_path.split("\\")[-3:-2])
+    path = "\\".join(s3_image_path.split("\\")[:-1])
+    date = " ".join(s3_image_path.split("\\")[-2:-1])
     date_s3 = date_s3 + date + "\n"
     
+    s3_image_path = s3_image_path + "\\xfdumanifest.xml"
+
     # Cut
     script_text = script_text + "echo \"\t Cut the image S3 " + date + "\"\n"
     cut_image = path + "\\S3_" + date + "_cut.dim"
-    text_cut = graph_cut.replace("!INPUT_Sentinel-3_LST!", s3_image_path).replace("!INPUT_AOI_WKT!", AOI_WTK).replace("!OUTPUT_Sentinel-3_cut!", cut_image)
-    path_cut = path_graph + "sentinel_3_cut_"  + date +".xml"
+    time_nc = path + "\\S3_" + date + "_time.nc"
+    text_cut = graph_cut.replace("!INPUT_Sentinel-3_LST!", s3_image_path).replace("!INPUT_AOI_WKT!", AOI_WTK).replace("!OUTPUT_Sentinel-3_cut!", cut_image).replace("!OUTPUT_Sentinel-3_time!", time_nc)
+    path_cut = path_graph + "sentinel_3_cut_time_"  + date +".xml"
     script_text = script_text + "gpt " + path_cut + "\n"
     f = open(path_cut, "w")
     f.write(text_cut)
@@ -109,10 +112,12 @@
 date_s2 = ""
 for line in file_path_s2:
     s2_image_path = line.strip()
-    path = "\\".join(s2_image_path.split("\\")[:-2])
-    date = " ".join(s2_image_path.split("\\")[-3:-2])
+    path = "\\".join(s2_image_path.split("\\")[:-1])
+    date = " ".join(s2_image_path.split("\\")[-2:-1])
     date_s2 = date_s2 + date + "\n"
     
+    s2_image_path = s2_image_path + "\\MTD_MSIL2A.xml"
+
     # Mask file
     script_text = script_text + "echo \"\t Computing the mask for the image S2 " + date + "\"\n"
     mask_image = path + "\\S2_" + date + "_mask.dim"
diff --git a/py_scripts/2_read_times_s3.py b/py_scripts/2_read_times_s3.py
@@ -0,0 +1,44 @@
+import xarray as xr
+import json
+
+path_param = "input\\parameters.json"
+path_file_date_s3 = "output\\days_s3.txt"
+path_file_datetime_s3 = "output\\days_time_s3.txt"
+
+# Definition of general varibles
+with open(path_param, "r") as f:
+    p = json.load(f)
+general_path = p["general_path"]
+
+# Variable for Sentinel 3 images
+file_date_s3 = open(path_file_date_s3, "r")
+date_s3 = []
+for line in file_date_s3:
+    date = line.strip()
+    date_s3.append(date)
+file_date_s3.close()
+
+text = ""
+for day in date_s3:
+    path = general_path + "S3\\"  + day + "\\S3_" + day + "_time.nc"
+    ds = xr.open_dataset(path)
+    start_time = int(ds.attrs["start_date"][12:14])*60 + int(ds.attrs["start_date"][15:17])
+    end_time = int(ds.attrs["stop_date"][12:14])*60 + int(ds.attrs["stop_date"][15:17])
+    mean_time = int((start_time + end_time)/2)
+    mean_time_h = int(mean_time/60)
+    mean_time_m = int(mean_time%60)
+    if mean_time_h < 10:
+        hours = "0" + str(mean_time_h)
+    else:
+        hours = str(mean_time_h)
+    if mean_time_m < 10:
+        minutes = "0" + str(mean_time_m)
+    else:
+        minutes = str(mean_time_m)
+    text = text + day + " " + hours + ":" + minutes + "\n"
+
+f = open(path_file_datetime_s3, "w")
+f.write(text)
+f.close()
+
+print("\tTimes file created!\n\t" + path_file_datetime_s3)
diff --git a/py_scripts/3_create_script_python.py b/py_scripts/3_create_script_python.py
@@ -29,6 +29,7 @@
 general_path_s3 = general_path + "S3\\YYYY_MM_DD\\S3_YYYY_MM_DD"
 general_path_era = general_path + "era\\YYYY_MM_DD\\S3_YYYY_MM_DD"
 general_path_out = general_path + "out\\YYYY_MM_DD\\out_YYYY_MM_DD"
+general_path_et = general_path + "et\\YYYY_MM_DD\\out_YYYY_MM_DD"
 general_path_s2 = general_path + "S2\\YYYY_MM_DD\\S2_YYYY_MM_DD"
 
 # Variable for Sentinel 2 images
@@ -61,7 +62,7 @@
 cmd_long = python_path + "\\python " + python_path + "\\sen-et-snap-scripts\\longwave_irradiance.py --meteo_product PATHERA_meteo.dim --at_band air_temperature --vp_band vapour_pressure --ap_band air_pressure --at_height 100 --output_file PATHOUT_longwave_ir.dim"
 cmd_short = python_path + "\\python " + python_path + "\\sen-et-snap-scripts\\net_shortwave_radiation.py --lsp_product PATHS2_leaf_spectra.dim --lai_product PATHS2_biophysical.dim --csp_product PATHS2_structural_parameters.dim --mi_product PATHERA_meteo.dim --sza_product PATHS3_hr_vza.dim --soil_ref_vis " + soil_ref_vis + " --soil_ref_nir " + soil_ref_nir + " --output_file PATHOUT_shortwave_ra.dim"
 cmd_fluxes = python_path + "\\python " + python_path + "\\sen-et-snap-scripts\\energy_fluxes.py --lst PATHS3_lst_sharpened.dim --lst_vza PATHS3_hr_vza.dim --lai PATHS2_biophysical.dim --csp PATHS2_structural_parameters.dim --fgv PATHS2_frac_green.dim --ar PATHS2_aerodynamic_parameters.dim --mi PATHERA_meteo.dim --nsr PATHOUT_shortwave_ra.dim --li PATHOUT_longwave_ir.dim --mask PATHS2_mask.dim --soil_roughness " + soil_roughness_en + " --alpha_pt " + alpha_pt + " --atmospheric_measurement_height 100 --green_vegetation_emissivity " + green_vegetation_emissivity + " --soil_emissivity " + soil_emissivity + " --save_component_fluxes 1 --save_component_temperature 1 --save_aerodynamic_parameters 1 --output_file PATHOUT_instantaneous_fluxes.dim"
-cmd_et = python_path + "\\python " + python_path + "\\sen-et-snap-scripts\\daily_evapotranspiration.py --ief_file PATHOUT_instantaneous_fluxes.dim --mi_file PATHERA_meteo.dim --output_file PATHOUT_daily_evapotranspiration.dim"
+cmd_et = python_path + "\\python " + python_path + "\\sen-et-snap-scripts\\daily_evapotranspiration.py --ief_file PATHOUT_instantaneous_fluxes.dim --mi_file PATHERA_meteo.dim --output_file PATHET_daily_evapotranspiration.dim"
 
 # creation of commands for S2 images
 text = ""
@@ -93,6 +94,7 @@
     t_general_path_s3 = general_path_s3.replace("YYYY", year).replace("MM", mon).replace("DD", day)
     t_general_path_era = general_path_era.replace("YYYY", year).replace("MM", mon).replace("DD", day)
     t_general_path_out = general_path_out.replace("YYYY", year).replace("MM", mon).replace("DD", day)
+    t_general_path_et = general_path_et.replace("YYYY", year).replace("MM", mon).replace("DD", day)
     data_test = datetime(int(year), int(mon), int(day))
     cloz_dict = {abs(data_test.timestamp() - date.timestamp()) : date for date in date_s2}
     close_s2 = cloz_dict[min(cloz_dict.keys())]
@@ -113,7 +115,7 @@
     text = text + "echo \"\t Computing the energy fluxes for the image S3 " + str(date) + "\"\n"
     text = text + cmd_fluxes.replace("PATHS3", t_general_path_s3).replace("PATHS2", t_general_path_s2).replace("PATHERA", t_general_path_era).replace("PATHOUT", t_general_path_out) + "\n"
     text = text + "echo \"\t Computing the evapotranspiration for the image S3 " + str(date) + "\"\n"
-    text = text + cmd_et.replace("PATHERA", t_general_path_era).replace("PATHOUT", t_general_path_out) + "\n"
+    text = text + cmd_et.replace("PATHERA", t_general_path_era).replace("PATHOUT", t_general_path_out).replace("PATHET", t_general_path_et) + "\n"
     text = text + "echo \"\t Finish the computation of the evapotranspiration for the image S3 " + str(date) + "\"\n\n"
 
 text = text.replace("\\", "\\\\")

Original file line number	Diff line number	Diff line change
`@@ -1 +1 @@`
`1`		`-D:\TEMP\S2\2019_08_06\S2A_MSIL2A_20190806T095031_N0213_R079_T33TWF_20190806T114240.SAFE\MTD_MSIL2A.xml`
	`1`	`+D:\TEMP\S2\2019_08_06\S2A_MSIL2A_20190806T095031_N0213_R079_T33TWF_20190806T114240.SAFE`
Original file line number	Diff line number	Diff line change
`@@ -1 +1 @@`
`1`		`-D:\TEMP\S3\2019_08_03\S3A_SL_2_LST____20190803T083233_20190803T101332_20190804T153249_6059_047_335______LN2_O_NT_003.SEN3\xfdumanifest.xml`
	`1`	`+D:\TEMP\S3\2019_08_03\S3A_SL_2_LST____20190803T083233_20190803T101332_20190804T153249_6059_047_335______LN2_O_NT_003.SEN3`