Add VIIRS obs processing to dev/gcafs.v1

dev/jobs/JGCDAS_PREPARE_OBS

@@ -0,0 +1,43 @@
+#! /usr/bin/env bash
+
+source "${HOMEgcafs}/ush/preamble.sh"
+source "${HOMEgcafs}/ush/jjob_header.sh"
+
+##############################################
+# Set variables used in the script
+##############################################
+
+##############################################
+# Begin JOB SPECIFIC work
+##############################################
+
+###############################################################
+# Run relevant script
+
+EXSCRIPT=${DUMPAODPY:-${HOMEgcafs}/scripts/exgcdas_prepare_obs.py}
+${EXSCRIPT}
+status=$?
+if [[ ${status} -ne 0 ]]; then
+    exit "${status}"
+fi
+
+##############################################
+# End JOB SPECIFIC work
+##############################################
+
+##############################################
+# Final processing
+##############################################
+if [[ -e "${pgmout}" ]]; then
+    cat "${pgmout}"
+fi
+
+##########################################
+# Remove the Temporary working directory
+##########################################
+cd "${DATAROOT}" || exit
+if [[ "${KEEPDATA}" == "NO" ]]; then
+    rm -rf "${DATA}"
+fi
+
+exit 0

dev/scripts/exgcdas_prepare_obs.py

@@ -0,0 +1,27 @@
+#!/usr/bin/env python3
+# exgcdas_prepare_obs.py
+# This script will collect and preprocess
+# aerosol optical depth observations for
+# global aerosol assimilation
+import os
+
+from wxflow import AttrDict, Logger, cast_strdict_as_dtypedict, parse_j2yaml
+from pyobsforge.task.aero_prepobs import AerosolObsPrep
+
+# Initialize root logger
+logger = Logger(level='DEBUG', colored_log=True)
+
+
+if __name__ == '__main__':
+
+    # Take configuration from environment and cast it as python dictionary
+    config_env = cast_strdict_as_dtypedict(os.environ)
+    # Take configuration from YAML file to augment/append config dict
+    config_yaml = parse_j2yaml(os.path.join(config_env['HOMEgcafs'], 'parm', 'chem', 'prepare_obs.yaml'), config_env)
+    # Combine configs together
+    config = AttrDict(**config_env, **config_yaml['aoddump'])
+
+    aeroObs = AerosolObsPrep(config)
+    aeroObs.initialize()
+    aeroObs.execute()
+    aeroObs.finalize()

dev/ush/setup_gcafs_for_nco.py

@@ -278,7 +278,7 @@ def setup_gcafs_for_nco():
     for file_path in all_copied_files:
         num_replacements = replace_gfs_with_gcafs(file_path)
         print(f"Modified {file_path}: {num_replacements} replacements made.")
-    
-    
+
+
 if __name__ == "__main__":
     setup_gcafs_for_nco()

parm/chem/prepare_obs.yaml

@@ -0,0 +1,9 @@
+aoddump:
+  provider: VIIRSAOD
+  platforms: ['npp', 'j01', 'n21'] # note j01==n20
+  binning_stride: 33
+  binning_min_number_of_obs: 32
+  binning_cressman_radius: 25
+  thinning_threshold: 0
+  channel: 4
+  preqc: 0

ush/python/pyobsforge/__init__.py

@@ -0,0 +1,5 @@
+import os
+
+__docformat__ = "restructuredtext"
+__version__ = "0.1.0"
+pyobsforge_directory = os.path.dirname(__file__)

ush/python/pyobsforge/obsdb/__init__.py

@@ -0,0 +1 @@
+from .obsdb import BaseDatabase  # noqa

ush/python/pyobsforge/obsdb/jrr_aod_db.py

@@ -0,0 +1,93 @@
+import os
+import glob
+from datetime import datetime, timedelta
+from pyobsforge.obsdb import BaseDatabase
+
+
+class JrrAodDatabase(BaseDatabase):
+    """Class to manage an observation file database for JRR-AOD data."""
+
+    def __init__(self, db_name="jrr_aod_obs.db",
+                 dcom_dir="/lfs/h1/ops/prod/dcom/",
+                 obs_dir="jrr_aod"):
+        base_dir = os.path.join(dcom_dir, '*', obs_dir)
+        super().__init__(db_name, base_dir)
+
+    def create_database(self):
+        """
+        Create the SQLite database and observation files table.
+
+        The table `obs_files` contains the following columns:
+        - `id`: A unique identifier for each record (auto-incremented primary key).
+        - `filename`: The full path to the observation file (must be unique).
+        - `obs_time`: The timestamp of the observation, extracted from the filename.
+        - `receipt_time`: The timestamp when the file was added to the `dcom` directory.
+        - `satellite`: The satellite from which the observation was collected (e.g., NPP, NOAA-20).
+        """
+        query = """
+        CREATE TABLE IF NOT EXISTS obs_files (
+            id INTEGER PRIMARY KEY AUTOINCREMENT,
+            filename TEXT UNIQUE,
+            obs_time TIMESTAMP,
+            receipt_time TIMESTAMP,
+            satellite TEXT
+        )
+        """
+        self.execute_query(query)
+
+    def parse_filename(self, filename):
+        """Extract metadata from filenames matching the JRR-AOD pattern."""
+        # Make sure the filename matches the expected pattern
+        # Pattern: JRR-AOD_v3r2_n21_sYYYYMMDDHHMMSS_eYYYYMMDDHHMMSS_cYYYYMMDDHHMMSS.nc
+        basename = os.path.basename(filename)
+        parts = basename.split('_')
+        try:
+            if len(parts) >= 4 and parts[0] == "JRR-AOD":
+                obs_time = datetime.strptime(parts[3][1:13], "%Y%m%d%H%M")
+                receipt_time = datetime.fromtimestamp(os.path.getctime(filename))
+                satellite = parts[2]
+                return filename, obs_time, receipt_time, satellite
+        except ValueError:
+            return None
+
+        return None
+
+    def ingest_files(self):
+        """Scan the directory for new JRR-AOD observation files and insert them into the database."""
+        obs_files = glob.glob(os.path.join(self.base_dir, "*.nc"))
+        print(f"Found {len(obs_files)} new files to ingest")
+
+        records_to_insert = []
+        for file in obs_files:
+            parsed_data = self.parse_filename(file)
+            if parsed_data:
+                records_to_insert.append(parsed_data)
+
+        if records_to_insert:
+            query = """
+                INSERT INTO obs_files (filename, obs_time, receipt_time, satellite)
+                VALUES (?, ?, ?, ?)
+            """
+            self.insert_records(query, records_to_insert)
+
+
+if __name__ == "__main__":
+    db = JrrAodDatabase(dcom_dir="/home/gvernier/Volumes/hera-s1/runs/realtimeobs/lfs/h1/ops/prod/dcom/")
+
+    # Check for new files
+    db.ingest_files()
+
+    # Query files for a given DA cycle
+    da_cycle = "20250316120000"
+    window_begin = datetime.strptime(da_cycle, "%Y%m%d%H%M%S") - timedelta(hours=3)
+    window_end = datetime.strptime(da_cycle, "%Y%m%d%H%M%S") + timedelta(hours=3)
+
+    valid_files = db.get_valid_files(window_begin=window_begin,
+                                     window_end=window_end)
+
+    print(f"Found {len(valid_files)} valid files for DA cycle {da_cycle}")
+    for valid_file in valid_files:
+        if os.path.exists(valid_file):
+            print(f"Valid file: {valid_file}")
+        else:
+            print(f"File does not exist: {valid_file}")

ush/python/pyobsforge/obsdb/obsdb.py

@@ -0,0 +1,144 @@
+from logging import getLogger
+import sqlite3
+from datetime import datetime, timedelta
+from wxflow.sqlitedb import SQLiteDB
+from wxflow import FileHandler
+from os.path import basename, join
+
+logger = getLogger(__name__.split('.')[-1])
+
+
+class BaseDatabase(SQLiteDB):
+    """Base class for managing different types of file-based databases."""
+
+    def __init__(self, db_name: str, base_dir: str) -> None:
+        """
+        Initialize the database.
+
+        :param db_name: Name of the SQLite database.
+        :param base_dir: Directory containing observation files.
+        """
+        super().__init__(db_name)
+        self.base_dir = base_dir
+        self.create_database()
+
+    def create_database(self):
+        """Create the SQLite database. Must be implemented by subclasses."""
+        raise NotImplementedError("Subclasses must implement create_database method")
+
+    def get_connection(self):
+        """Return the database connection."""
+        return self.connection
+
+    def parse_filename(self):
+        """Parse a filename and extract relevant metadata. Must be implemented by subclasses."""
+        raise NotImplementedError("Subclasses must implement parse_filename method")
+
+    def ingest_files(self):
+        """Scan the directory for new observation files and insert them into the database."""
+        raise NotImplementedError("Subclasses must implement ingest_files method")
+
+    def insert_record(self, query: str, params: tuple) -> None:
+        """Insert a record into the database."""
+        self.connect()
+        cursor = self.connection.cursor()
+        try:
+            cursor.execute(query, params)
+            self.connection.commit()
+        except sqlite3.IntegrityError:
+            pass  # Skip duplicates
+        finally:
+            self.disconnect()
+
+    def insert_records(self, query: str, params_list: list[tuple]) -> None:
+        """
+        Insert multiple records into the database.
+
+        :param query: SQL query for inserting records.
+        :param params_list: List of tuples containing the parameters for each record.
+        """
+        self.connect()
+        cursor = self.connection.cursor()
+        try:
+            cursor.executemany(query, params_list)
+            self.connection.commit()
+        except sqlite3.IntegrityError:
+            pass  # Skip duplicates
+        finally:
+            self.disconnect()
+
+    def execute_query(self, query: str, params: tuple = None) -> list:
+        """Execute a query and return the results."""
+        self.connect()
+        cursor = self.connection.cursor()
+        cursor.execute(query, params or [])
+        results = cursor.fetchall()
+        self.disconnect()
+        return results
+
+    def get_valid_files(self,
+                        window_begin: datetime,
+                        window_end: datetime,
+                        dst_dir: str,
+                        instrument: str = None,
+                        satellite: str = None,
+                        obs_type: str = None,
+                        check_receipt: str = "none") -> list:
+        """
+        Retrieve and copy to dst_dir a list of observation files within a specified time window, possibly filtered by instrument,
+        satellite, and observation type. The check_receipt parameter can be 'gdas', 'gfs', or 'none'. If 'gdas' or
+        'gfs' is specified, files are further filtered based on their receipt time to ensure they meet the
+        required delay criteria.
+
+        :param window_begin: Start of the time window (datetime object).
+        :param window_end: End of the time window (datetime object).
+        :param dst_dir: Destination directory where valid files will be copied.
+        :param instrument: (Optional) Filter by instrument name.
+        :param satellite: (Optional) Filter by satellite name.
+        :param obs_type: (Optional) Filter by observation type.
+        :param check_receipt: (Optional) Specify receipt time check ('gdas', 'gfs', or 'none').
+        :return: List of valid observation file paths in the destination directory.
+        """
+
+        query = """
+        SELECT filename FROM obs_files
+        WHERE obs_time BETWEEN ? AND ?
+        """
+        minutes_behind_realtime = {'gdas': 160, 'gfs': 20}
+        params = [window_begin, window_end]
+
+        if instrument:
+            query += " AND instrument = ?"
+            params.append(instrument)
+        if satellite:
+            query += " AND satellite = ?"
+            params.append(satellite)
+        if obs_type:
+            query += " AND obs_type = ?"
+            params.append(obs_type)
+
+        results = self.execute_query(query, tuple(params))
+        valid_files = []
+        for row in results:
+            filename = row[0]
+            if check_receipt in ["gdas", "gfs"]:
+                query = "SELECT receipt_time FROM obs_files WHERE filename = ?"
+                receipt_time = self.execute_query(query, (filename,))[0][0]
+                receipt_time = datetime.strptime(receipt_time, "%Y-%m-%d %H:%M:%S.%f")
+                if receipt_time <= window_end - timedelta(minutes=minutes_behind_realtime[check_receipt]):
+                    continue
+
+            valid_files.append(filename)
+
+        # Copy files to the destination directory
+        dst_files = []
+        if len(valid_files) > 0:
+            src_dst_obs_list = []  # list of [src_file, dst_file]
+            for src_file in valid_files:
+                dst_file = join(dst_dir, f"{basename(src_file)}")
+                dst_files.append(dst_file)
+                src_dst_obs_list.append([src_file, dst_file])
+            FileHandler({'mkdir': [dst_dir]}).sync()
+            FileHandler({'copy': src_dst_obs_list}).sync()
+
+        return dst_files