From ee7862915cb02c1226ecb9a060daf52cfa98257f Mon Sep 17 00:00:00 2001
From: BaptisteVandecrux <b.vandecrux@gmail.com>
Date: Wed, 28 Sep 2022 13:56:11 +0200
Subject: [PATCH] now prints lat lon in fortran input file

---
 sice_f.py | 17 +++++++++++++----
 1 file changed, 13 insertions(+), 4 deletions(-)

diff --git a/sice_f.py b/sice_f.py
index 3cb5423..7dd5b44 100644
--- a/sice_f.py
+++ b/sice_f.py
@@ -17,6 +17,7 @@
 #import bav_lib as bl
 import pandas as pd
 import rioxarray
+import xarray as xr
 import shutil
 import glob
 import subprocess
@@ -27,6 +28,14 @@
 start_time = time.time()
 InputFolder = InputFolder + "/"
 Oa01 = rio.open(InputFolder + "r_TOA_01.tif")
+band1 = Oa01.read(1)
+height = band1.shape[0]
+width = band1.shape[1]
+cols, rows = np.meshgrid(np.arange(width), np.arange(height))
+xs, ys = rio.transform.xy(Oa01.transform, rows, cols)
+lons= np.array(xs)
+lats = np.array(ys)
+print('lons shape', lons.shape)
 meta = Oa01.meta
 
 def WriteOutput(var, var_name, in_folder):
@@ -56,10 +65,10 @@ def WriteOutput(var, var_name, in_folder):
 
 olci_toa = np.vstack(
     (
-        np.arange(1, len(sza.flatten()) + 1),  # pixel number_x
-        np.arange(1, len(sza.flatten()) + 1),  # pixel number_y
-        np.arange(1, len(sza.flatten()) + 1),  # latitude
-        np.arange(1, len(sza.flatten()) + 1),  # longitude
+        cols.flatten(),  # pixel number_x
+        rows.flatten(),  # pixel number_y
+        lats.flatten(),  # latitude
+        lons.flatten(),  # longitude
         sza.flatten(),  # solar zenith angle
         saa.flatten(),  # soalr azimuthal angle
         vza.flatten(),  # viewing zenith angle