Merge pull request #52 from C2SM-RCM/coli

Updates from new developpement

docs/source/api.rst

@@ -64,7 +64,12 @@ the anthropogenic and biogenic part of a CO2.
 
 .. autofunction:: emiproc.speciation.read_speciation_table
 
+.. autofunction:: emiproc.speciation.merge_substances
 
+Input/Output
+^^^^^^^^^^^^
+
+.. autofunction:: emiproc.speciation.read_speciation_table
 
 Utilities
 ---------

docs/source/installation.rst

@@ -25,22 +25,6 @@ the cloned code.
 
     pip install -e .
 
-Many packages will be installed as dependencies.
 
 
-On Windows 
-----------
 
-emiproc works also on windows, but many dependencies are not built for 
-windows.
-
-However you can use `pipwin` which will let you install all the dependencies
-for windows in the same way `pip` works.
-
-.. code::
-    
-    pip install pipwin 
-
-    pipwin install geos
-    pipwin install geopandas
-    pipwin install cartopy 

emiproc/exports/icon.py

@@ -366,7 +366,7 @@ def make_icon_time_profiles(
                         dt_end,
                         time_profiles[key],
                         local_tz=tz,
-                        freq="H",
+                        freq="h",
                     ).to_numpy()  # [max_shift + shift : -max_shift + shift - 1]
                     for shift, tz in zip(shifts, time_zones)
                 ]
@@ -378,7 +378,7 @@ def make_icon_time_profiles(
                 coords={
                     "datetime": (
                         ("hourofyear",),
-                        pd.date_range(dt_start, dt_end, freq="H"),
+                        pd.date_range(dt_start, dt_end, freq="h"),
                     ),
                     "timezone_of_country": (("country",), time_zones),
                 },

emiproc/exports/rasters.py

@@ -7,7 +7,7 @@
 from emiproc.grids import RegularGrid
 from emiproc.regrid import remap_inventory
 from emiproc.exports.netcdf import NetcdfAttributes
-from emiproc.utilities import Units, SEC_PER_YR
+from emiproc.utilities import Units, SEC_PER_YR, PER_CELL_UNITS, PER_M2_UNITS
 
 
 def export_raster_netcdf(
@@ -78,6 +78,10 @@ def export_raster_netcdf(
         ) * 1e9
     else:
         raise NotImplementedError(f"Unknown {unit=}")
+
+    unit_str = str(unit.value)
+    if unit in PER_CELL_UNITS:
+        unit_str += " cell-1"
 
     ds = xr.Dataset(
         data_vars=(
@@ -89,7 +93,7 @@ def export_raster_netcdf(
                     {
                         "standard_name": f"{sub}_{cat}",
                         "long_name": f"{sub}_{cat}",
-                        "units": str(unit.value),
+                        "units": unit_str,
                         "comment": f"emissions of {sub} in {cat}",
                         "projection": f"{crs}",
                     },
@@ -119,7 +123,7 @@ def export_raster_netcdf(
                     {
                         "standard_name": f"tendency_of_atmosphere_mass_content_of_{sub}_due_to_emission",
                         "long_name": f"{sub}",
-                        "units": str(unit.value),
+                        "units": unit_str,
                         "comment": f"emissions of {sub}",
                         "projection": f"{crs}",
                     },
@@ -139,7 +143,7 @@ def export_raster_netcdf(
                     "long_name": "longitude",
                     "units": "degrees_east",
                     "comment": "center_of_cell",
-                    "bounds": "lon_bnds",
+                    "bounds": f"{lon_name}_bnds",
                     "projection": f"{crs}",
                     "axis": "X",
                 },
@@ -152,11 +156,34 @@ def export_raster_netcdf(
                     "long_name": "latitude",
                     "units": "degrees_north",
                     "comment": "center_of_cell",
-                    "bounds": "lat_bnds",
+                    "bounds": f"{lat_name}_bnds",
                     "projection": f"{crs}",
                     "axis": "Y",
                 },
             ),
+            f"{lon_name}_bnds": (
+                f"{lon_name}_bnds",
+                grid.lon_bounds,
+                {
+                    "standard_name": "longitude",
+                    "long_name": "longitude of cell boundaries",
+                    "units": "degrees_east",
+                    "comment": "cell_boundaries",
+                    "axis": "X",
+                },
+            ),
+            f"{lat_name}_bnds": (
+                f"{lat_name}_bnds",
+                grid.lat_bounds,
+                {
+                    "standard_name": "latitude",
+                    "long_name": "latitude of cell boundaries",
+                    "units": "degrees_north",
+                    "comment": "cell_boundaries",
+                    "axis": "Y",
+                },
+            ),
+
         },
         attrs=netcdf_attributes,
     )
@@ -178,7 +205,7 @@ def export_raster_netcdf(
             ds[f"emi_{sub}_all_sectors"].attrs = {
                 "standard_name": f"tendency_of_atmosphere_mass_content_of_{sub}_due_to_emission",
                 "long_name": f"Aggregated Emissions of {sub} from all sectors",
-                "units": str(unit.value),
+                "units": unit_str,
                 "comment": "annual mean emission rate",
                 "projection": f"{crs}",
             }
@@ -218,19 +245,19 @@ def export_raster_netcdf(
                 "comment": "annual total emission",
             }
 
-    if unit in [Units.KG_PER_M2_PER_S, Units.MUG_PER_M2_PER_S]:
-        # add the cell area
-        ds["cell_area"] = (
-            [lat_name, lon_name],
-            np.array(grid.cell_areas).reshape(grid.shape).T,
-            {
-                "standard_name": "cell_area",
-                "long_name": "cell_area",
-                "units": "m2",
-                "comment": "area of the cell",
-                "projection": f"{crs}",
-            },
-        )
+
+    # add the cell area
+    ds["cell_area"] = (
+        [lat_name, lon_name],
+        np.array(grid.cell_areas).reshape(grid.shape).T,
+        {
+            "standard_name": "cell_area",
+            "long_name": "cell_area",
+            "units": "m2",
+            "comment": "area of the cell",
+            "projection": f"{crs}",
+        },
+    )
     path = Path(path)
     out_filepath = path.with_suffix(".nc")
     ds.to_netcdf(out_filepath)

emiproc/grids.py

@@ -51,6 +51,9 @@ def __init__(self, name: str, crs: int | str = WGS84):
         """
         self.name = name
         self.crs = crs
+
+    def __str__(self) -> str:
+        return f"{self.__class__.__name__}({self.name})"
 
     @property
     def gdf(self) -> gpd.GeoDataFrame:
@@ -158,6 +161,10 @@ class RegularGrid(Grid):
     lon_range: np.ndarray
     lat_range: np.ndarray
 
+    # The edges of the cells
+    lat_bounds: np.ndarray
+    lon_bounds: np.ndarray
+
     xmin: float
     xmax: float
     ymin: float
@@ -214,6 +221,7 @@ def __init__(
         if xmax is None and ymax is None:
             xmax = xmin + nx * dx
             ymax = ymin + ny * dy
+        self.xmax, self.ymax = xmax, ymax
 
         # Calculate all grid parameters
         self.nx, self.ny = nx, ny
@@ -222,8 +230,16 @@ def __init__(
         self.lon_range = np.arange(xmin, xmax, self.dx) + self.dx / 2
         self.lat_range = np.arange(ymin, ymax, self.dy) + self.dy / 2
 
+        self.lon_bounds = np.concatenate([self.lon_range - self.dx / 2, [self.lon_range[-1] + self.dx / 2]])
+        self.lat_bounds = np.concatenate([self.lat_range - self.dy / 2, [self.lat_range[-1] + self.dy / 2]])
+
+        assert len(self.lon_range) == nx
+        assert len(self.lat_range) == ny
+        assert len(self.lon_bounds) == nx + 1
+        assert len(self.lat_bounds) == ny + 1
+
         if name is None:
-            name = f"RegularGrid_x({xmin},{xmax})_y({ymin},{ymax})_nx({nx})_ny({ny})"
+            name = f"x({xmin},{xmax})_y({ymin},{ymax})_nx({nx})_ny({ny})"
 
         super().__init__(name, crs)
 

emiproc/inventories/__init__.py

@@ -239,7 +239,11 @@ def substances(self) -> list[Substance]:
 
     @property
     def total_emissions(self) -> pd.DataFrame:
-        """Simple accessor to the function."""
+        """Calculate the total emissions, returning a DataFrame.
+        
+        Simple accessor to the function
+        :py:func:`~emiproc.inventories.utils.get_total_emissions`.
+        """
         from emiproc.inventories.utils import get_total_emissions
 
         return pd.DataFrame(get_total_emissions(self)).T
@@ -421,6 +425,23 @@ def set_profile(
             if "category" not in indexes_array.dims:
                 # Exapnd the array over the categories of the inventory
                 indexes_array = indexes_array.expand_dims({"category": self.categories})
+            if category not in indexes_array.coords["category"].values:
+                if category not in self.categories:
+                    raise ValueError(
+                        f"Category {category} is not in the inventory. "
+                        "Please add it before setting the profile."
+                    )
+                # Otherwise add it to the coords and set the values to -1 (unassigned)
+                # emtpy da with the new category being the only one in the 
+                new_cat_array = xr.DataArray(
+                    -1,
+                    dims=indexes_array.dims,
+                    coords={
+                        "category": [category],
+                        **{coord: indexes_array.coords[coord] for coord in indexes_array.dims if coord != "category"},
+                    },
+                )
+                indexes_array = xr.concat([indexes_array, new_cat_array], dim="category")
 
             sel_dict["category"] = category
         if substance is not None:

emiproc/inventories/categories_groups.py

@@ -8,9 +8,9 @@
     ],
     # Industry
     "GNFR_B": [
-        "eipro",
+        "eipro",  # Flächenquellen Industrie
         "eipwp",  # this is the weitere punktquelle (additional point sources)
-        "eipzm",
+        "eipzm",  # Punktquellen Zementwerke (point sources cement plants)
     ],
     # Other stationary combustion (services, residential, agriculture)
     "GNFR_C": [
@@ -24,10 +24,10 @@
     ],
     # Fugitives
     "GNFR_D": [
-        "eilgk", # F-Gase: Läden, Gebäude mit Klimaanlagen Refrigeration (F-gases: shops, buildings with air conditioning)
-        "eivgn", # Verluste Gasnetz (losses gas network)
-        "evklm", # F-Gase: Klimaanlagen Motorfahrzeuge (F-gases: air conditioning motor vehicles)
-        "evtrk", # F-Gase: Transporte mit Kühlung  (F-gases: transports with cooling)
+        "eilgk",  # F-Gase: Läden, Gebäude mit Klimaanlagen Refrigeration (F-gases: shops, buildings with air conditioning)
+        "eivgn",  # Verluste Gasnetz (losses gas network)
+        "evklm",  # F-Gase: Klimaanlagen Motorfahrzeuge (F-gases: air conditioning motor vehicles)
+        "evtrk",  # F-Gase: Transporte mit Kühlung  (F-gases: transports with cooling)
     ],
     # Solvents and product use
     "GNFR_E": [
@@ -44,8 +44,8 @@
         # "evstr_nmvoc",
         # "evstr_nox",
         # "evstr_so2",
-        "evstr", # Strassenverkehr (road transport)
-        "evzon", # Zonenverkehr Kaltstart/Verdampfung (zone traffic, cold start/evaporation)
+        "evstr",  # Strassenverkehr (road transport)
+        "evzon",  # Zonenverkehr Kaltstart/Verdampfung (zone traffic, cold start/evaporation)
     ],
     # Shipping
     "GNFR_G": [

emiproc/inventories/utils.py

@@ -679,6 +679,7 @@ def drop(
         inv: Inventory,
         substances: list[Substance] = [],
         categories: list[Category] = [],
+        keep_instead_of_drop: bool = False,
 ) -> Inventory:
     """Drop substances and categories from an inventory.
 
@@ -691,6 +692,8 @@ def drop(
     :arg inv: The inventory to drop from.
     :arg substances: The substances to drop.
     :arg categories: The categories to drop.
+    :arg keep_instead_of_drop: If True, the substances and categories will be kept
+        instead of being dropped.
     """
 
     # Check the types 
@@ -699,6 +702,13 @@ def drop(
             raise TypeError(f"{var_name=} should be a list.")
         if not all(isinstance(sub, str) for sub in var):
             raise TypeError(f"{var_name=} should be a list of strings.")
+
+    if keep_instead_of_drop:
+        # Keep instead of drop
+        if substances:
+            substances = [sub for sub in inv.substances if sub not in substances]
+        if categories:
+            categories = [cat for cat in inv.categories if cat not in categories]
 
     # Deep copy of the inventory
     out_inv = inv.copy(no_gdfs=True)