Merge branch 'devel' of gitlab.in2p3.fr:satinv/cif into devel

pycif/plugins/models/chimere/io/inputs/make_obs.py

@@ -57,8 +57,10 @@ def make_obs(self, datastore, runsubdir, mode, tracer, do_simu=True):
             ddata = d[1]
             spec = ddata["parameter"].upper()
 
-            # convert level to fortran
+            # convert level and coordinates to fortran
             ddata["level"] += 1
+            ddata["i"] += 1
+            ddata["j"] += 1
 
             # TODO: deal with altitude
             if np.isnan(ddata["level"]) or ddata["level"] <= 0:

pycif/plugins/obsvects/standard/utils/hcoord.py

@@ -43,6 +43,7 @@ def hcoord(obsvect, chunksize=1e3, **kwargs):
     is_unstructured = getattr(obsvect.domain, "unstructured_domain",
                               False)
     
+    chunksize = 2e6
     # Makes simplified operations if regular
     if isregular or is_unstructured:
         lon = ds["lon"]
@@ -114,7 +115,7 @@ def hcoord(obsvect, chunksize=1e3, **kwargs):
     # Cropping observations outside the model domain
     mask = ~np.isnan(ds["i"]) & ~np.isnan(ds["j"])
     obsvect.datastore["data"] = ds.loc[mask]
-
+    
     return obsvect
 
 
@@ -155,6 +156,10 @@ def find_gridcell(
                 zlonc[0, np.newaxis] - lon[:, np.newaxis], discont=discont
             )
             xlat = 1.0 / (lat[:, np.newaxis] - zlatc[:, 0])
+            
+            xlon[xlon == np.inf] = -np.inf
+            xlat[xlat == np.inf] = -np.inf
+            
             i = xlat.argmin(axis=1)
             j = xlon.argmin(axis=1)
             
@@ -164,6 +169,14 @@ def find_gridcell(
             maskj = (lon < zlonc.min()) | (lon > zlonc.max())
             i = np.where(maski | maskj, np.nan, i)
             j = np.where(maski | maskj, np.nan, j)
+            
+            # Points in last column or row should be shifted
+            # as corners were used
+            maski = (i >= nlat)
+            maskj = (i >= nlon)
+            i = np.where(maski, nlat - 1, i)
+            j = np.where(maskj, nlon - 1, j)
+            
             return i, j
 
         # If only one (lon, lat) tuple was provided

pycif/plugins/transforms/satellites/forward.py

@@ -44,11 +44,14 @@ def forward(
     ddi = min(di, df)
     ref_indexes = ~target_datastore.duplicated(subset=["indorig"])
     y0 = target_datastore.loc[ref_indexes]
-
+    
+    # Exit if empty observation datastore
+    if target_datastore.size == 0:
+        return
+    
     # Building the extended dataframe
     iq1 = y0["station"]
     list_satIDs = iq1.unique()
-
     ds_p = target_datastore.set_index("indorig")[
         ["pressure", "dp", "airm", "hlay"]]
     for satID in list_satIDs:

pycif/plugins/transforms/satellites/vinterp.py

@@ -44,7 +44,7 @@ def vertical_interp(pres_in, pres_out, cropstrato):
     levout = scipy.interpolate.griddata(points_in, levin, points_out).reshape(
         (-1, lenout)
     )
-
+    
     # Getting interpolation coefficients
     xlow = np.floor(levout).astype(int)
     xhigh = xlow + 1

pycif/plugins/transforms/toobsvect/adjoint.py

@@ -97,7 +97,8 @@ def adjoint(
             obsvect.dy[param.ypointer: param.ypointer + param.dim][mask]
         datastore[tracer_id]["data"].loc[:, "sim"] = \
             obsvect.ysim[param.ypointer: param.ypointer + param.dim][mask]
-        datastore[tracer_id]["data"].loc[:, "sim_tl"] = 0.
+        datastore[tracer_id]["data"].loc[:, "sim_tl"] = \
+            0. * obsvect.ysim[param.ypointer: param.ypointer + param.dim][mask]
         
         # Initializing the local domain for this set of observations
         if not getattr(param, "obs_domain_init", False):