Merge branch 'main-dev' into investigate1291

README.rst

@@ -24,5 +24,5 @@ pyaerocom
 .. |Coverage| image:: https://codecov.io/gh/metno/pyaerocom/branch/main-dev/graph/badge.svg?token=A0AdX8YciZ
     :target: https://codecov.io/gh/metno/pyaerocom
 
-.. |DOI| image:: https://zenodo.org/badge/DOI/10.5281/zenodo.10374181.svg
-   :target: https://doi.org/10.5281/zenodo.10374181
+.. |DOI| image:: https://zenodo.org/badge/DOI/10.5281/zenodo.13270713.svg
+  :target: https://doi.org/10.5281/zenodo.13270713

docs/_static/aeroval/cfg_examples_example1.py

@@ -179,7 +179,7 @@
 that this does note cover all possible settings, but only the ones that
 deviate from the default co-location settings.
 For all co-location options (and their defaults), see here:
-https://pyaerocom.readthedocs.io/en/latest/api.html?highlight=ColocationSetup#pyaerocom.colocation_auto.ColocationSetup
+https://pyaerocom.readthedocs.io/en/latest/api.html#pyaerocom.colocation.colocation_setup.ColocationSetup
 
 Note that all of the available co-location settings can also be specified
 for model and obs entries individually below.

pyaerocom/aeroval/_processing_base.py

@@ -158,7 +158,7 @@ class DataImporter(HasColocator):
     are / can be specified flexibly for each model and obs entry in an
     analysis setup (:class:`EvalSetup`). Proper handling of these reading
     constraints and data import settings are handled in the
-    :class:`pyaerocom.colocation_auto.Colocator` engine, therefore the reading
+    :class:`pyaerocom.colocation.Colocator` engine, therefore the reading
     in this class is done via the :class:`Colocator` engine.
 
 

pyaerocom/aeroval/data/var_scale_colmap.ini

@@ -438,3 +438,6 @@ colmap = coolwarm
 scale = [0.0, 0.004, 0.008, 0.012, 0.016, 0.02, 0.04, 0.06, 0.08, 0.1, 0.2, 0.3, 0.4]
 colmap = coolwarm
 
+[vcdno2]
+scale = [0, 10, 20, 30, 40, 50, 60, 70, 80, 90, 100, 110, 120]
+colmap = coolwarm

pyaerocom/aeroval/data/var_web_info.ini

@@ -802,3 +802,23 @@ menu_name = ratio PM10 PM2.5
 vertical_type = 3D
 category = Particle ratio
 
+# Vertical column densities
+[vcdno2]
+menu_name = VCD NO2
+vertical_type = 3D
+category = Vertical column density
+
+[vcdso2]
+menu_name = VCD SO2
+vertical_type = 3D
+category = Vertical column density
+
+[vcdhcho]
+menu_name = VCD HCHO
+vertical_type = 3D
+category = Vertical column density
+
+[vcdco]
+menu_name = VCD CO
+vertical_type = 3D
+category = Vertical column density

pyaerocom/aeroval/experiment_processor.py

@@ -24,7 +24,7 @@ class ExperimentProcessor(ProcessingEngine, HasColocator):
     processing engine will perform spatial and temporal co-location and will
     store on co-located NetCDF file (e.g. if there are 2 models, 2 observation
     networks and 2 variables there will be 4 co-located NetCDF files).
-    The co-location is done using :class:`pyaerocom.colocation_auto.Colocator`.
+    The co-location is done using :class:`pyaerocom.Colocator`.
 
     """
 
@@ -73,7 +73,9 @@ def _run_single_entry(self, model_name, obs_name, var_list):
                 var_list = list(set(obs_vars) & set(var_list))
                 if not var_list:
                     logger.warning(
-                        "var_list %s and obs_vars %s mismatch.", var_list_asked, obs_vars
+                        "var_list %s and obs_vars %s mismatch.",
+                        var_list_asked,
+                        obs_vars,
                     )
                     return
 

pyaerocom/aeroval/glob_defaults.py

@@ -133,6 +133,7 @@ class CategoryType(str, Enum):
     deposition = "Deposition"
     temperature = "Temperature"
     particle_ratio = "Particle ratio"
+    vertical_column_density = "Vertical column density"
     UNDEFINED = "UNDEFINED"
 
     def __str__(self):

pyaerocom/aeroval/setupclasses.py

@@ -138,14 +138,14 @@ class StatisticsSetup(BaseModel, extra="allow"):
     annual_stats_constrained : bool
         if True, then only sites are considered that satisfy a potentially
         specified annual resampling constraint (see
-        :attr:`pyaerocom.colocation_auto.ColocationSetup.min_num_obs`). E.g.
+        :attr:`pyaerocom.colocation.ColocationSetup.min_num_obs`). E.g.
 
         lets say you want to calculate statistics (bias,
         correlation, etc.) for monthly model / obs data for a given site and
         year. Lets further say, that there are only 8 valid months of data, and
         4 months are missing, so statistics will be calculated for that year
         based on 8 vs. 8 values. Now if
-        :attr:`pyaerocom.colocation_auto.ColocationSetup.min_num_obs` is
+        :attr:`pyaerocom.colocation.ColocationSetup.min_num_obs` is
         specified in way that requires e.g. at least 9 valid months to
         represent the whole year, then this station will not be considered in
         case `annual_stats_constrained` is True, else it will. Defaults to
@@ -457,7 +457,11 @@ def colocation_opts(self) -> ColocationSetup:
             if key in ColocationSetup.model_fields
         }
         # need to pass some default values to the ColocationSetup if not provided in config
-        default_dict = {"save_coldata": True, "keep_data": False, "resample_how": "mean"}
+        default_dict = {
+            "save_coldata": True,
+            "keep_data": False,
+            "resample_how": "mean",
+        }
         for key in default_dict:
             if key not in model_args:
                 model_args[key] = default_dict[key]

pyaerocom/colocation/colocation_utils.py

@@ -394,7 +394,14 @@ def check_ts_type(data, ts_type):
 
 
 def _colocate_site_data_helper(
-    stat_data, stat_data_ref, var, var_ref, ts_type, resample_how, min_num_obs, use_climatology_ref
+    stat_data,
+    stat_data_ref,
+    var,
+    var_ref,
+    ts_type,
+    resample_how,
+    min_num_obs,
+    use_climatology_ref,
 ):
     """
     Helper method that colocates two timeseries from 2 StationData objects
@@ -445,15 +452,26 @@ def _colocate_site_data_helper(
         obs_ts = stat_data_ref.calc_climatology(var_ref, min_num_obs=min_num_obs)[var_ref]
     else:
         obs_ts = stat_data_ref.resample_time(
-            var_ref, ts_type=ts_type, how=resample_how, min_num_obs=min_num_obs, inplace=True
+            var_ref,
+            ts_type=ts_type,
+            how=resample_how,
+            min_num_obs=min_num_obs,
+            inplace=True,
         )[var_ref]
 
     # fill up missing time stamps
     return pd.concat([obs_ts, grid_ts], axis=1, keys=["ref", "data"])
 
 
 def _colocate_site_data_helper_timecol(
-    stat_data, stat_data_ref, var, var_ref, ts_type, resample_how, min_num_obs, use_climatology_ref
+    stat_data,
+    stat_data_ref,
+    var,
+    var_ref,
+    ts_type,
+    resample_how,
+    min_num_obs,
+    use_climatology_ref,
 ):
     """
     Helper method that colocates two timeseries from 2 StationData objects
@@ -514,7 +532,11 @@ def _colocate_site_data_helper_timecol(
     # =============================================================================
 
     stat_data.resample_time(
-        var_name=var, ts_type=str(coltst), how=resample_how, min_num_obs=min_num_obs, inplace=True
+        var_name=var,
+        ts_type=str(coltst),
+        how=resample_how,
+        min_num_obs=min_num_obs,
+        inplace=True,
     )
 
     stat_data_ref.resample_time(
@@ -594,8 +616,8 @@ def colocate_gridded_ungridded(
 ):
     """Colocate gridded with ungridded data (low level method)
 
-    For high-level colocation see :class:`pyaerocom.colocation_auto.Colocator`
-    and :class:`pyaerocom.colocation_auto.ColocationSetup`
+    For high-level colocation see :class:`pyaerocom.colocation.Colocator`
+    and :class:`pyaerocom.ColocationSetup`
 
     Note
     ----
@@ -996,7 +1018,10 @@ def correct_model_stp_coldata(coldata, p0=None, t0=273.15, inplace=False):
     arr = coldata.data
 
     coords = zip(
-        arr.latitude.values, arr.longitude.values, arr.altitude.values, arr.station_name.values
+        arr.latitude.values,
+        arr.longitude.values,
+        arr.altitude.values,
+        arr.station_name.values,
     )
     if p0 is None:
         p0 = pressure()  # STD conditions sea level

pyaerocom/config.py

@@ -102,6 +102,9 @@ class Config:
     #: ICOS name
     ICOS_NAME = "ICOS"
 
+    # TROPOMI access names
+    TROPOMI_XEMEP_R01x01_NAME = "TROPOMI_XEMEP_R01x01"
+
     #: boolean specifying wheter EBAS DB is copied to local cache for faster
     #: access, defaults to True
     EBAS_DB_LOCAL_CACHE = True
@@ -185,7 +188,11 @@ class Config:
 
     # this dictionary links environment ID's with corresponding subdirectory
     # names that are required to exist in order to load this environment
-    _check_subdirs_cfg = {"metno": "aerocom", "users-db": "AMAP", "local-db": "modeldata"}
+    _check_subdirs_cfg = {
+        "metno": "aerocom",
+        "users-db": "AMAP",
+        "local-db": "modeldata",
+    }
 
     with resources.path("pyaerocom.data", "variables.ini") as path:
         _var_info_file = str(path)
@@ -751,7 +758,11 @@ def reload(self, keep_basedirs=True):
         self.read_config(self.last_config_file, keep_basedirs)
 
     def read_config(
-        self, config_file, basedir=None, init_obslocs_ungridded=False, init_data_search_dirs=False
+        self,
+        config_file,
+        basedir=None,
+        init_obslocs_ungridded=False,
+        init_data_search_dirs=False,
     ):
         """
         Import paths from one of the config ini files

pyaerocom/data/paths.ini

@@ -102,6 +102,8 @@ ICOS = ${BASEDIR}/aerocom/aerocom1/AEROCOM_OBSDATA/ICOS/aggregated/
 #ICPFORESTS = /lustre/storeB/project/fou/kl/emep/People/danielh/projects/pyaerocom/obs/ICP-forests/dep/
 ICPFORESTS = ${BASEDIR}/aerocom/aerocom1/AEROCOM_OBSDATA/icp-forests/dep/
 
+TROPOMI_XEMEP_R01x01=/lustre/storeB/project/fou/kl/sesam/work/CSO-gridded/xEMEP__r01x01/data/
+
 [obsnames]
 #names of the different obs networks
 #Aeronet V3
@@ -137,6 +139,7 @@ CNEMC = CNEMC
 ICOS = ICOS
 
 ICPFORESTS = ICPFORESTS
+TROPOMI_XEMEP_R01x01 = TROPOMI_XEMEP_R01x01
 
 [parameters]
 #parameters definition
@@ -174,4 +177,5 @@ EARLINET = 2000
 EEA_NRT = 2020
 EEA_V2 = 2016
 MEP = 2013
-ICOS = 2008
+ICOS = 2008
+TROPOMI_XEMEP_R01x01 = 2018

pyaerocom/data/variables.ini

@@ -5968,9 +5968,46 @@ unit = ug m-3
 description=Mass concentration of coarse elemental carbon
 unit = ug C m-3x
 
-
-
-
 [conclevoglucosan]
 description=Mass concentration of Levoglucosan
 unit = ug m-3
+
+[vcdno2]
+var_name = vcdno2
+description = Vertical Column Density NO2
+standard_name = troposphere_mole_content_of_nitrogen_dioxide
+var_type = vertical column density
+unit = umol m-2
+minimum = 0
+maximum = 100000
+dimensions = time,lat,lon
+
+[vcdso2]
+var_name = vcdso2
+description = Vertical Column Density SO2
+standard_name = atmosphere_mole_content_of_sulfur_dioxide
+var_type = vertical column density
+unit = mol m-2
+minimum = 0
+maximum = 100000000
+dimensions = time,lat,lon
+
+[vcdhcho]
+var_name = vcdhcho
+description = Vertical Column Density HCHO
+standard_name = troposphere_mole_content_of_formaldehyde
+var_type = vertical column density
+unit = mol m-2
+minimum = 0
+maximum = 100000000
+dimensions = time,lat,lon
+
+[vcdco]
+var_name = vcdco
+description = Vertical Column Density CO
+standard_name = atmosphere_mole_content_of_carbon_monoxide
+var_type = vertical column density
+unit = mol m-2
+minimum = 0
+maximum = 100000000
+dimensions = time,lat,lon

pyaerocom/helpers.py

@@ -945,14 +945,14 @@ def merge_station_data(
     return merged
 
 
-def _get_pandas_freq_and_loffset(freq):
+def _get_pandas_freq_and_offset(freq: str) -> tuple[str, pd.Timedelta | None]:
     """Helper to convert resampling info"""
     if freq in TS_TYPE_TO_PANDAS_FREQ:
         freq = TS_TYPE_TO_PANDAS_FREQ[freq]
-    loffset = None
+    offset = None
     if freq in PANDAS_RESAMPLE_OFFSETS:
-        loffset = PANDAS_RESAMPLE_OFFSETS[freq]
-    return (freq, loffset)
+        offset = PANDAS_RESAMPLE_OFFSETS[freq]
+    return (freq, offset)
 
 
 def make_datetime_index(start, stop, freq):
@@ -983,10 +983,10 @@ def make_datetime_index(start, stop, freq):
     if not isinstance(stop, pd.Timestamp):
         stop = to_pandas_timestamp(stop)
 
-    freq, loffset = _get_pandas_freq_and_loffset(freq)
+    freq, offset = _get_pandas_freq_and_offset(freq)
     idx = pd.date_range(start=start, end=stop, freq=freq)
-    if loffset is not None:
-        idx = idx + pd.Timedelta(loffset)
+    if offset is not None:
+        idx = idx + offset
     return idx
 
 
@@ -1096,7 +1096,7 @@ def resample_timeseries(ts, freq, how=None, min_num_obs=None):
         p = int(how.split("percentile")[0])
         how = lambda x: np.nanpercentile(x, p)  # noqa: E731
 
-    freq, loffset = _get_pandas_freq_and_loffset(freq)
+    freq, offset = _get_pandas_freq_and_offset(freq)
     resampler = ts.resample(freq)
 
     data = resampler.agg(how)
@@ -1106,8 +1106,8 @@ def resample_timeseries(ts, freq, how=None, min_num_obs=None):
         invalid = numobs < min_num_obs
         if np.any(invalid):
             data.values[invalid] = np.nan
-    if loffset is not None:
-        data.index = data.index + pd.Timedelta(loffset)
+    if offset is not None:
+        data.index = data.index + offset
     return data
 
 
@@ -1165,8 +1165,8 @@ def resample_time_dataarray(arr, freq, how=None, min_num_obs=None):
     if min_num_obs is not None:
         invalid = arr.resample(time=pd_freq).count(dim="time") < min_num_obs
 
-    freq, loffset = _get_pandas_freq_and_loffset(freq)
-    resampler = arr.resample(time=pd_freq, loffset=loffset)
+    freq, offset = _get_pandas_freq_and_offset(freq)
+    resampler = arr.resample(time=pd_freq, offset=offset)
     try:
         aggfun = getattr(resampler, how)
     except AttributeError:

pyaerocom/stats/mda8/mda8.py

@@ -104,9 +104,6 @@ def _rolling_average_8hr(arr: xr.DataArray) -> xr.DataArray:
 
 
 def _daily_max(arr: xr.DataArray) -> xr.DataArray:
-    # TODO: Base is deprecated, and using offset="1h" is the proper way to do this.
-    # However this currently breaks the old-dependencies test in CI. Should be
-    # changed in the future.
-    return arr.resample(time="24H", base=1).reduce(
+    return arr.resample(time="24H", offset="1h").reduce(
         lambda x, axis: np.apply_along_axis(min_periods_max, 1, x, min_periods=18)
     )