NOA-ReACT / PollyXT-SCC-Pipelines

Tool for the processing of PollyXT using Single Calculus Chain
https://noa-react.github.io/PollyXT-SCC-Pipelines/
GNU Lesser General Public License v3.0
3 stars 2 forks source link

depol channel @1064nm #26

Closed ulysses78 closed 1 year ago

ulysses78 commented 1 year ago

Hey @thgeorgiou , at Mindelo (Cape Verde) we have a PollyXT system running, which has a 355, 532 and also a 1064nm depol-channel. Locally at our server I managed to add the missing lines to the pollyxt_pipelines script. The following scripts have to be changed:

locations/__init__.py
polly_to_scc/scc_netcdf.py

Adding the following parameters:

total_channel_1064_nm_idx
cross_channel_1064_nm_idx
calibration_1064nm_total_channel_ids
calibration_1064nm_cross_channel_ids

Of course one has to change every occurency and initialization of it. I will add my code down here.

One thing to keep in mind in the locations.ini file: at the moment the scripts are written in a way, that one has to add for every site a "pseudo" 1064nm depol channel, which is kind of ugly. Is there a way to avoid this? Could you please implement this missing 1064nm depol-channel to your official scripts and create a pip-release?

FYI: @HolgerPollyNet

ulysses78 commented 1 year ago

Here is the code for the locations/__init__.py:

"""
Contains information about locations

Each location (i.e. SCC station) is defined in an .ini file. For some stations,
the .ini files are included with the software but custom locations can be defined.
"""

import io
import sys
from configparser import ConfigParser, SectionProxy
from importlib.resources import read_text
from typing import Dict, List, NamedTuple, Union

from rich.markdown import Markdown
from rich.table import Table

from pollyxt_pipelines import config
from pollyxt_pipelines.console import console
from pollyxt_pipelines.utils import ints_to_csv

class Location(NamedTuple):
    """
    Represents a physical location of PollyXT installation.
    """

    name: str
    """Location friendly name"""

    profile_name: str
    """How are the WRF profile names prefixed"""

    sounding_provider: str
    """Which radiosonde provider to use"""

    scc_code: str
    """SCC Station code"""

    lat: float
    """Latitude of station"""

    lon: float
    """Longitude of station"""

    altitude_asl: float
    """Altitude of station"""

    daytime_configuration: int
    """SCC Lidar Configuration ID - Daytime"""

    nighttime_configuration: int
    """SCC Lidar Configuration ID - Nightime"""

    calibration_configuration_355nm: int
    """SCC Lidar Configuration ID - Calibration (355 nm)"""

    calibration_configuration_532nm: int
    """SCC Lidar Configuration ID - Calibration (532 nm)"""

    calibration_configuration_1064nm: int
    """SCC Lidar Configuration ID - Calibration (1064 nm)"""

    depol_calibration_zero_state: int
    """Value of `depol_cal_angle` when there is *no* calibration taking place"""

    channel_id: List[int]
    """Mapping of PollyXT Channels to SCC Channels
    Comma-separated list. The order of the list is the order of the channels in the
    PollyXT netCDF file.
    """

    background_low: List[int]
    """Value for the `Background_Low` variable"""

    background_high: List[int]
    """Value for the `Background_High` variable"""

    lr_input: List[int]
    """Value for the `lr_input` variable"""

    temperature: int
    """Temperature at the lidar station (`Temperature_at_Lidar_Station` variable)"""

    pressure: int
    """Pressure at the lidar station (`Pressure_at_Lidar_Station` variable)"""

    total_channel_355_nm_idx: int
    """Index in Polly netCDF file for the total channel (355nm)"""

    cross_channel_355_nm_idx: int
    """Index in Polly netCDF file for the cross channel (355nm)"""

    total_channel_532_nm_idx: int
    """Index in Polly netCDF file for the total channel (532nm)"""

    cross_channel_532_nm_idx: int
    """Index in Polly netCDF file for the cross channel (532nm)"""

    total_channel_1064_nm_idx: int
    """Index in Polly netCDF file for the total channel (1064nm)"""

    cross_channel_1064_nm_idx: int
    """Index in Polly netCDF file for the cross channel (1064nm)"""

    calibration_355nm_total_channel_ids: List[int]
    """
    Calibration channel SCC IDs for 355nm. Comma separated list. First value must be the
    +45° channel, second value must be the -45° channel.
    """

    calibration_355nm_cross_channel_ids: List[int]
    """
    Calibration channel SCC IDs for 355nm. Comma separated list. First value must be the
    +45° channel, second value must be the -45° channel.
    """

    calibration_532nm_total_channel_ids: List[int]
    """
    Calibration channel SCC IDs for 532nm. Comma separated list. First value must be the
    +45° channel, second value must be the -45° channel.
    """

    calibration_532nm_cross_channel_ids: List[int]
    """
    Calibration channel SCC IDs for 532nm. Comma separated list. First value must be the
    +45° channel, second value must be the -45° channel.
    """

    calibration_1064nm_total_channel_ids: List[int]
    """
    Calibration channel SCC IDs for 1064nm. Comma separated list. First value must be the
    +45° channel, second value must be the -45° channel.
    """

    calibration_1064nm_cross_channel_ids: List[int]
    """
    Calibration channel SCC IDs for 1064nm. Comma separated list. First value must be the
    +45° channel, second value must be the -45° channel.
    """

    def print(self):
        """
        Prints this location as a Table in the terminal
        """

        table = Table(title=self.name)
        table.add_column("Key")
        table.add_column("Value")

        for key, value in self._asdict().items():
            if isinstance(value, list):
                value = ints_to_csv(value)
            table.add_row(key, str(value))

        console.print(table)

def location_from_section(name: str, section: SectionProxy) -> Location:
    """
    Create a Location from a ConfigParser Section (SectionProxy)
    """

    channel_id = [int(x.strip()) for x in section.get("channel_id").split(",")]
    background_low = [int(x.strip()) for x in section.get("background_low").split(",")]
    background_high = [
        int(x.strip()) for x in section.get("background_high").split(",")
    ]
    lr_input = [int(x.strip()) for x in section.get("lr_input").split(",")]

    calibration_355nm_total_channel_ids = [
        int(x.strip())
        for x in section.get("calibration_355nm_total_channel_ids").split(",")
    ]
    calibration_355nm_cross_channel_ids = [
        int(x.strip())
        for x in section.get("calibration_355nm_cross_channel_ids").split(",")
    ]
    calibration_532nm_total_channel_ids = [
        int(x.strip())
        for x in section.get("calibration_532nm_total_channel_ids").split(",")
    ]
    calibration_532nm_cross_channel_ids = [
        int(x.strip())
        for x in section.get("calibration_532nm_cross_channel_ids").split(",")
    ]
    calibration_1064nm_total_channel_ids = [
        int(x.strip())
        for x in section.get("calibration_1064nm_total_channel_ids").split(",")
    ]
    calibration_1064nm_cross_channel_ids = [
        int(x.strip())
        for x in section.get("calibration_1064nm_cross_channel_ids").split(",")
    ]

    return Location(
        name=name,
        scc_code=section["scc_code"],
        lat=section.getfloat("lat"),
        lon=section.getfloat("lon"),
        altitude_asl=section.getfloat("altitude_asl"),
        daytime_configuration=section.getint("daytime_configuration"),
        nighttime_configuration=section.getint("nighttime_configuration"),
        calibration_configuration_355nm=section.getint(
            "calibration_configuration_355nm"
        ),
        calibration_configuration_532nm=section.getint(
            "calibration_configuration_532nm"
        ),
        calibration_configuration_1064nm=section.getint(
            "calibration_configuration_1064nm"
        ),
        depol_calibration_zero_state=section.getint(
            "depol_calibration_zero_state",
        ),
        channel_id=channel_id,
        background_low=background_low,
        background_high=background_high,
        lr_input=lr_input,
        temperature=section.getint("temperature"),
        pressure=section.getint("pressure"),
        total_channel_355_nm_idx=section.getint("total_channel_355_nm_idx"),
        cross_channel_355_nm_idx=section.getint("cross_channel_355_nm_idx"),
        total_channel_532_nm_idx=section.getint("total_channel_532_nm_idx"),
        cross_channel_532_nm_idx=section.getint("cross_channel_532_nm_idx"),
        total_channel_1064_nm_idx=section.getint("total_channel_1064_nm_idx"),
        cross_channel_1064_nm_idx=section.getint("cross_channel_1064_nm_idx"),
        calibration_355nm_total_channel_ids=calibration_355nm_total_channel_ids,
        calibration_355nm_cross_channel_ids=calibration_355nm_cross_channel_ids,
        calibration_532nm_total_channel_ids=calibration_532nm_total_channel_ids,
        calibration_532nm_cross_channel_ids=calibration_532nm_cross_channel_ids,
        calibration_1064nm_total_channel_ids=calibration_1064nm_total_channel_ids,
        calibration_1064nm_cross_channel_ids=calibration_1064nm_cross_channel_ids,
        sounding_provider=section["sounding_provider"],
        profile_name=section["profile_name"],
    )

def read_locations() -> Dict[str, Location]:
    """
    Reads all built-in and custom locations into a dictionary: name -> Location
    """

    locations = {}

    # Read built-in locations
    locations_buffer = io.StringIO(
        read_text("pollyxt_pipelines.locations", "locations.ini")
    )
    locations_config = ConfigParser()
    locations_config.read_file(locations_buffer)

    for name in locations_config.sections():
        section = locations_config[name]
        locations[name] = location_from_section(name, section)

    # Read custom locations
    location_paths = [path / "locations.ini" for path in config.config_paths()]
    locations_config = ConfigParser()
    locations_config.read(location_paths)

    for name in locations_config.sections():
        section = locations_config[name]
        try:
            locations[name] = location_from_section(name, section)
        except Exception:
            console.print(
                f"Could not load locations from config file, problem occured in section [{name}]."
            )
            console.print("Check the following files:")
            for path in location_paths:
                if path.is_file():
                    console.print(f"\t-{path}")
            sys.exit(1)

    return locations

LOCATIONS = read_locations()
"""List of all known locations"""

def get_location_by_scc_code(code: str) -> Union[Location, None]:
    """
    Returns a location by its SCC code or `None` if it doesn't exist.
    """

    for loc in LOCATIONS.values():
        if loc.scc_code == code:
            return loc
    return None

def unknown_location_error(name: str):
    """
    Prints an error message that the given location is not found, along with a
    list of known locations
    """
    error = (
        f"[error]Could not find location[/error]{name}[error]\nKnown locations:\n\n."
    )
    for l in LOCATIONS:
        error += f"* {l.name}"

    console.print(Markdown(error))
ulysses78 commented 1 year ago

Here is the code for the polly_to_scc/scc_netcdf.py:


"""
Routines for converting PollyXT files to SCC files
"""

from datetime import timedelta
from enum import Enum, IntEnum
from pathlib import Path
from typing import Tuple

import numpy as np
from netCDF4 import Dataset

from pollyxt_pipelines import utils
from pollyxt_pipelines.locations import Location
from pollyxt_pipelines.polly_to_scc import pollyxt
from pollyxt_pipelines.polly_to_scc.exceptions import (
    NoMeasurementsInTimePeriod,
    TimeOutsideFile,
)

class Wavelength(Enum):
    """Laser wavelength"""

    NM_355 = 355
    NM_532 = 532
    NM_1064 = 1064

class Atmosphere(IntEnum):
    """
    Which atmosphere to use to processing. These are the possible values for `molecular_calc`.
    """

    AUTOMATIC = 0
    RADIOSONDE = 1
    CLOUDNET = 2
    STANDARD_ATMOSPHERE = 4

    @staticmethod
    def from_string(x: str):
        x = x.lower().strip()
        if x == "automatic":
            return Atmosphere.AUTOMATIC
        if x == "radiosonde":
            return Atmosphere.RADIOSONDE
        if x == "cloudnet":
            return Atmosphere.CLOUDNET
        if x == "standard":
            return Atmosphere.STANDARD_ATMOSPHERE

        raise ValueError(f"Unknown atmosphere {x}")

def create_scc_netcdf(
    pf: pollyxt.PollyXTFile,
    output_path: Path,
    location: Location,
    atmosphere=Atmosphere.STANDARD_ATMOSPHERE,
) -> Tuple[str, Path]:
    """
    Convert a PollyXT netCDF file to a SCC file.

    Parameters:
        pf: An opened PollyXT file. When you create this, you can specify the time period of interest.
        output_path: Where to store the produced netCDF file
        location: Where did this measurement take place
        atmosphere: What kind of atmosphere to use.

    Note:
        If atmosphere is set to Atmosphere.SOUNDING, the `Sounding_File_Name` attribute will be set to
        `rs_{MEASUREMENT_ID}.rs`, ie the filename of the accompaning radiosonde. This file is *not* created
        by this function.

    Returns:
        A tuple containing  the measurement ID and the output path
    """

    # Calculate measurement ID
    measurement_id = pf.start_date.strftime(f"%Y%m%d{location.scc_code}%H%M")

    # Create SCC file
    # Output filename is always the measurement ID
    output_filename = output_path / f"{measurement_id}.nc"
    nc = Dataset(output_filename, "w")

    # Create dimensions (mandatory!)
    nc.createDimension("points", np.size(pf.raw_signal, axis=1))
    nc.createDimension("channels", np.size(pf.raw_signal, axis=2))
    nc.createDimension("time", None)
    nc.createDimension("nb_of_time_scales", 1)
    nc.createDimension("scan_angles", 1)

    # Create Global Attributes (mandatory!)
    nc.Measurement_ID = measurement_id
    nc.RawData_Start_Date = pf.start_date.strftime("%Y%m%d")
    nc.RawData_Start_Time_UT = pf.start_date.strftime("%H%M%S")
    nc.RawData_Stop_Time_UT = pf.end_date.strftime("%H%M%S")

    # Create Global Attributes (optional)
    nc.RawBck_Start_Date = nc.RawData_Start_Date
    nc.RawBck_Start_Time_UT = nc.RawData_Start_Time_UT
    nc.RawBck_Stop_Time_UT = nc.RawData_Stop_Time_UT
    if atmosphere == Atmosphere.RADIOSONDE:
        nc.Sounding_File_Name = f"rs_{measurement_id[:-2]}.nc"
    # nc.Overlap_File_Name = 'ov_' + selected_start.strftime('%Y%m%daky%H') + '.nc'

    # Custom attribute for configuration ID
    # From 04:00 until 16:00 we use daytime configuration
    if pf.start_date.replace(
        hour=4, minute=0
    ) < pf.start_date and pf.start_date < pf.start_date.replace(hour=16, minute=0):
        nc.X_PollyXTPipelines_Configuration_ID = location.daytime_configuration
    else:
        nc.X_PollyXTPipelines_Configuration_ID = location.nighttime_configuration

    # Create Variables. (mandatory)
    raw_data_start_time = nc.createVariable(
        "Raw_Data_Start_Time", "i4", dimensions=("time", "nb_of_time_scales"), zlib=True
    )
    raw_data_stop_time = nc.createVariable(
        "Raw_Data_Stop_Time", "i4", dimensions=("time", "nb_of_time_scales"), zlib=True
    )
    raw_lidar_data = nc.createVariable(
        "Raw_Lidar_Data", "f8", dimensions=("time", "channels", "points"), zlib=True
    )
    channel_id = nc.createVariable(
        "channel_ID", "i4", dimensions=("channels"), zlib=True
    )
    id_timescale = nc.createVariable(
        "id_timescale", "i4", dimensions=("channels"), zlib=True
    )
    laser_pointing_angle = nc.createVariable(
        "Laser_Pointing_Angle", "f8", dimensions=("scan_angles"), zlib=True
    )
    laser_pointing_angle_of_profiles = nc.createVariable(
        "Laser_Pointing_Angle_of_Profiles",
        "i4",
        dimensions=("time", "nb_of_time_scales"),
        zlib=True,
    )
    laser_shots = nc.createVariable(
        "Laser_Shots", "i4", dimensions=("time", "channels"), zlib=True
    )
    background_low = nc.createVariable(
        "Background_Low", "f8", dimensions=("channels"), zlib=True
    )
    background_high = nc.createVariable(
        "Background_High", "f8", dimensions=("channels"), zlib=True
    )
    molecular_calc = nc.createVariable("Molecular_Calc", "i4", dimensions=(), zlib=True)
    nc.createVariable("Pol_Calib_Range_Min", "f8", dimensions=("channels"), zlib=True)
    nc.createVariable("Pol_Calib_Range_Max", "f8", dimensions=("channels"), zlib=True)
    pressure_at_lidar_station = nc.createVariable(
        "Pressure_at_Lidar_Station", "f8", dimensions=(), zlib=True
    )
    temperature_at_lidar_station = nc.createVariable(
        "Temperature_at_Lidar_Station", "f8", dimensions=(), zlib=True
    )
    lr_input = nc.createVariable("LR_Input", "i4", dimensions=("channels"), zlib=True)

    # Fill Variables with Data. (mandatory)
    raw_data_start_time[:] = (
        pf.measurement_time[~pf.calibration_mask, 1] - pf.measurement_time[0, 1]
    )
    raw_data_stop_time[:] = (
        pf.measurement_time[~pf.calibration_mask, 1] - pf.measurement_time[0, 1]
    ) + 30
    raw_lidar_data[:] = pf.raw_signal_swap[~pf.calibration_mask]
    channel_id[:] = np.array(location.channel_id)
    id_timescale[:] = np.zeros(np.size(pf.raw_signal[~pf.calibration_mask], axis=2))
    laser_pointing_angle[:] = int(pf.zenith_angle.item(0))
    laser_pointing_angle_of_profiles[:] = np.zeros(
        np.size(pf.raw_signal[~pf.calibration_mask], axis=0)
    )
    laser_shots[:] = pf.measurement_shots[~pf.calibration_mask]
    background_low[:] = np.array(location.background_low)
    background_high[:] = np.array(location.background_high)
    molecular_calc[:] = int(atmosphere)
    pressure_at_lidar_station[:] = location.pressure
    temperature_at_lidar_station[:] = location.temperature
    lr_input[:] = np.array(location.lr_input)

    # Close the netCDF file.
    nc.close()

    return measurement_id, output_filename

def create_scc_calibration_netcdf(
    pf: pollyxt.PollyXTFile,
    output_path: Path,
    location: Location,
    wavelength: Wavelength,
    pol_calib_range_min: int = 1200,
    pol_calib_range_max: int = 2500,
) -> Tuple[str, Path]:
    """
    From a PollyXT netCDF file, create the corresponding calibration SCC file.
    Calibration only occures when `depol_cal_angle` is not equal to the default state value.
    Take care to create the `PollyXTFile` with these intervals.

    Parameters:
        pf: An opened PollyXT file
        output_path: Where to store the produced netCDF file
        location: Where did this measurement take place
        wavelength: Calibration for 355nm or 532nm
        pol_calib_range_min: Calibration contant calculation, minimum height
        pol_calib_range_max: Calibration contant calculation, maximum height

    Returns:
        A tuple containing the measurement ID and the output path
    """

    # Calculate measurement ID
    measurement_id = pf.start_date.strftime(f"%Y%m%d{location.scc_code}%H")

    # Create SCC file
    # Output filename is always the measurement ID
    if wavelength == Wavelength.NM_355:
        output_filename = output_path / f"calibration_{measurement_id}_355.nc"
    elif wavelength == Wavelength.NM_532:
        output_filename = output_path / f"calibration_{measurement_id}_532.nc"
    elif wavelength == Wavelength.NM_1064:
        output_filename = output_path / f"calibration_{measurement_id}_1064.nc"
    else:
        raise ValueError(f"Unknown wavelength {wavelength}")
    nc = Dataset(output_filename, "w")

    # Find start/end indices for the +45 and -45 degree calibration cycles in Polly file
    idx = list(np.where(np.diff(pf.depol_cal_angle))[0])
    start_positive = 2
    idx = list(filter(lambda x: x >= start_positive + 4, idx))
    end_positive = idx[0]
    positive_length = end_positive - start_positive

    start_negative = idx[0] + 3
    idx = list(filter(lambda x: x >= start_negative + 4, idx))
    end_negative = pf.depol_cal_angle.shape[0] - 3
    negative_length = end_negative - start_negative

    # Reduce the larger period to match
    if positive_length > negative_length:
        end_positive -= positive_length - negative_length
        positive_length = negative_length
    elif negative_length > positive_length:
        end_negative -= negative_length - positive_length
        negative_length = positive_length

    # Create Dimensions. (mandatory)
    nc.createDimension("points", np.size(pf.raw_signal, axis=1))
    nc.createDimension("channels", 4)
    nc.createDimension("time", positive_length)
    nc.createDimension("nb_of_time_scales", 1)
    nc.createDimension("scan_angles", 1)

    # Create Global Attributes. (mandatory)
    # Move start date a couple of profiles forward to accomodate the fact that we skip
    # some profiles at the beginning of the file.
    start_date = pf.start_date + timedelta(seconds=(start_positive * 30))
    nc.RawData_Start_Date = start_date.strftime("%Y%m%d")
    nc.RawData_Start_Time_UT = start_date.strftime("%H%M%S")
    nc.RawData_Stop_Time_UT = pf.end_date.strftime("%H%M%S")

    # Create Global Attributes (optional)
    nc.RawBck_Start_Date = nc.RawData_Start_Date
    nc.RawBck_Start_Time_UT = nc.RawData_Start_Time_UT
    nc.RawBck_Stop_Time_UT = nc.RawData_Stop_Time_UT

    # Create Variables. (mandatory)
    raw_data_start_time = nc.createVariable(
        "Raw_Data_Start_Time", "i4", dimensions=("time", "nb_of_time_scales"), zlib=True
    )
    raw_data_stop_time = nc.createVariable(
        "Raw_Data_Stop_Time", "i4", dimensions=("time", "nb_of_time_scales"), zlib=True
    )
    raw_lidar_data = nc.createVariable(
        "Raw_Lidar_Data", "f8", dimensions=("time", "channels", "points"), zlib=True
    )
    channel_id = nc.createVariable(
        "channel_ID", "i4", dimensions=("channels"), zlib=True
    )
    id_timescale = nc.createVariable(
        "id_timescale", "i4", dimensions=("channels"), zlib=True
    )
    laser_pointing_angle = nc.createVariable(
        "Laser_Pointing_Angle", "f8", dimensions=("scan_angles"), zlib=True
    )
    laser_pointing_angle_of_profiles = nc.createVariable(
        "Laser_Pointing_Angle_of_Profiles",
        "i4",
        dimensions=("time", "nb_of_time_scales"),
        zlib=True,
    )
    laser_shots = nc.createVariable(
        "Laser_Shots", "i4", dimensions=("time", "channels"), zlib=True
    )
    background_low = nc.createVariable(
        "Background_Low", "f8", dimensions=("channels"), zlib=True
    )
    background_high = nc.createVariable(
        "Background_High", "f8", dimensions=("channels"), zlib=True
    )
    molecular_calc = nc.createVariable("Molecular_Calc", "i4", dimensions=(), zlib=True)
    pol_calib_range_min_var = nc.createVariable(
        "Pol_Calib_Range_Min", "f8", dimensions=("channels"), zlib=True
    )
    pol_calib_range_max_var = nc.createVariable(
        "Pol_Calib_Range_Max", "f8", dimensions=("channels"), zlib=True
    )
    pressure_at_lidar_station = nc.createVariable(
        "Pressure_at_Lidar_Station", "f8", dimensions=(), zlib=True
    )
    temperature_at_lidar_station = nc.createVariable(
        "Temperature_at_Lidar_Station", "f8", dimensions=(), zlib=True
    )

    # Fill Variables with Data. (mandatory)
    raw_data_start_time[:] = (
        pf.measurement_time[start_positive:end_positive, 1]
        - pf.measurement_time[start_positive, 1]
    )
    raw_data_stop_time[:] = (
        pf.measurement_time[start_negative:end_negative, 1]
        - pf.measurement_time[start_positive, 1]
    )
    id_timescale[:] = np.array([0, 0, 0, 0])
    laser_pointing_angle[:] = 5
    laser_pointing_angle_of_profiles[:, :] = 0.0
    laser_shots[:] = 600
    background_low[:] = np.array([0, 0, 0, 0])
    background_high[:] = np.array([249, 249, 249, 249])
    molecular_calc[:] = 0
    pol_calib_range_min_var[:] = np.repeat(pol_calib_range_min, 4)
    pol_calib_range_max_var[:] = np.repeat(pol_calib_range_max, 4)
    pressure_at_lidar_station[:] = location.pressure
    temperature_at_lidar_station[:] = location.temperature

    # Define total and cross channels IDs from Polly
    if wavelength == Wavelength.NM_355:
        total_channel_idx = location.total_channel_355_nm_idx
        cross_channel_idx = location.cross_channel_355_nm_idx
        channel_id[:] = np.array(
            location.calibration_355nm_total_channel_ids
            + location.calibration_355nm_cross_channel_ids
        )
        nc.Measurement_ID = measurement_id + "35"
        nc.X_PollyXTPipelines_Configuration_ID = (
            location.calibration_configuration_355nm
        )
    elif wavelength == Wavelength.NM_532:
        total_channel_idx = location.total_channel_532_nm_idx
        cross_channel_idx = location.cross_channel_532_nm_idx
        channel_id[:] = np.array(
            location.calibration_532nm_total_channel_ids
            + location.calibration_532nm_cross_channel_ids
        )
        nc.Measurement_ID = measurement_id + "53"
        nc.X_PollyXTPipelines_Configuration_ID = (
            location.calibration_configuration_532nm
        )
    elif wavelength == Wavelength.NM_1064:
        total_channel_idx = location.total_channel_1064_nm_idx
        cross_channel_idx = location.cross_channel_1064_nm_idx
        channel_id[:] = np.array(
            location.calibration_1064nm_total_channel_ids
            + location.calibration_1064nm_cross_channel_ids
        )
        nc.Measurement_ID = measurement_id + "10"
        nc.X_PollyXTPipelines_Configuration_ID = (
            location.calibration_configuration_1064nm
        )
    else:
        raise ValueError(f"Unknown wavelength {wavelength}")

    raw_lidar_data[:] = 0
    # Total channel, +45°
    raw_lidar_data[:, 0, :] = pf.raw_signal_swap[
        start_positive:end_positive, total_channel_idx, :
    ]
    # Cross channel, +45°
    raw_lidar_data[:, 2, :] = pf.raw_signal_swap[
        start_positive:end_positive, cross_channel_idx, :
    ]
    # Total channel, -45°
    raw_lidar_data[:, 1, :] = pf.raw_signal_swap[
        start_negative:end_negative, total_channel_idx, :
    ]
    # Cross channel, -45°
    raw_lidar_data[:, 3, :] = pf.raw_signal_swap[
        start_negative:end_negative, cross_channel_idx, :
    ]

    # Close the netCDF file.
    nc.close()

    return measurement_id, output_filename

def convert_pollyxt_file(
    repo: pollyxt.PollyXTRepository,
    output_path: Path,
    location: Location,
    interval: timedelta,
    atmosphere: Atmosphere,
    should_round=False,
    calibration=True,
    start_time=None,
    end_time=None,
):
    """
    Converts a pollyXT repository into a collection of SCC files. The input files will be split/merged into intervals
    before being converted to the new format.

    This function is a generator, so you can use it in a for loop to monitor progress:

        for measurement_id, path, start_time, end_time in convert_pollyxt_file(...):
            # Do something with id/path, maybe print a message?

    Parameters:
        repo: PollyXT file to convert
        output_path: Directory to write the SCC files
        location: Geographical information, where the measurement took place
        interval: What interval to use when splitting the PollyXT file (e.g. 1 hour)
        atmosphere: Which atmosphere to use on SCC
        should_round: If true, the interval starts will be rounded down. For example, from 01:02 to 01:00.
        calibration: Set to False to disable generation of calibration files.
        start_hour: Optionally, set when the first file should start. The intervals will start from here. (HH:MM or YYYY-MM-DD_HH:MM format, string)
        end_hour: Optionally, also set the end time. Must be used with `start_hour`. If this is set, only one output file
                  is generated, for your target interval (HH:MM or YYYY-MM-DD_HH:MM format, string).
    """

    # Open input netCDF
    measurement_start, measurement_end = repo.get_time_period()

    # Handle start/end time
    if start_time is not None:
        start_time = utils.date_option_to_datetime(measurement_start, start_time)

        if start_time < measurement_start or measurement_end < start_time:
            raise TimeOutsideFile(measurement_start, measurement_end, start_time)

        measurement_start = start_time

    if start_time is None and end_time is not None:
        raise ValueError("Can't use end_hour without start_hour")

    if end_time is not None:
        end_time = utils.date_option_to_datetime(measurement_end, end_time)

        if end_time < measurement_start or measurement_end < start_time:
            raise TimeOutsideFile(measurement_start, measurement_end, end_time)

        measurement_end = end_time
        interval = timedelta(seconds=(end_time - start_time).total_seconds())

    # Create output files
    interval_start = measurement_start
    while interval_start < measurement_end:
        # If the option is set, round down hours
        if should_round:
            interval_start = interval_start.replace(microsecond=0, second=0, minute=0)

        # Interval end
        interval_end = interval_start + interval

        # Open netCDF file and convert to SCC
        try:
            pf = repo.get_pollyxt_file(
                interval_start, interval_end + timedelta(seconds=30)
            )
            id, path = create_scc_netcdf(pf, output_path, location, atmosphere)

            yield id, path, pf.start_date, pf.end_date
        except NoMeasurementsInTimePeriod as ex:
            # Go to next loop
            interval_start = interval_end
            continue

        # Set start of next interval to the end of this one
        interval_start = interval_end

    # Generate calibration files
    if calibration:
        # Check for any valid calibration intervals
        for start, end in repo.get_calibration_periods():
            if start > measurement_start and end < measurement_end:
                pf = repo.get_pollyxt_file(start, end)

                id, path = create_scc_calibration_netcdf(
                    pf, output_path, location, wavelength=Wavelength.NM_532
                )
                yield id, path, start, end

                id, path = create_scc_calibration_netcdf(
                    pf, output_path, location, wavelength=Wavelength.NM_355
                )
                yield id, path, start, end

                id, path = create_scc_calibration_netcdf(
                    pf, output_path, location, wavelength=Wavelength.NM_1064
                )
                yield id, path, start, end
ulysses78 commented 1 year ago

Just search within these two scripts for "1064" and you will get all new entries ;-)

elmarinou commented 1 year ago

hello Holger! we have a new colleague Evgenia in cc. She is taking over the polly automations from Thanasis (guided by him), as he is shifted 100% to l2a+ assimilation. so you can mention her for the optimizations discussions!

Evgenia, take action!!

Holger, looking foreword for suggestions for upgrades!

Best wishes,

Eleni

Στις 2023-02-07 09:05, ulysses78 έγραψε:

Hey @thgeorgiou [1] , at Mindelo (Cape Verde) we have a PollyXT system running, which has a 355, 532 and also a 1064nm depol-channel. Locally at our server I managed to add the missing lines to the pollyxt_pipelines script. The following scripts have to be changed:

locations/init.py polly_to_scc/scc_netcdf.py

Adding the following parameters:

total_channel_1064_nm_idx cross_channel_1064_nm_idx calibration_1064nm_total_channel_ids calibration_1064nm_cross_channel_ids

Of course one has to change every occurency and initialization of it. I will add my code down here.

One thing to keep in mind in the locations.ini file: at the moment the scripts are written in a way, that one has to add for every site a "pseudo" 1064nm depol channel, which is kind of ugly. Is there a way to avoid this? Could you please implement this missing 1064nm depol-channel to your official scripts and create a pip-release?

FYI: @HolgerPollyNet [2]

-- Reply to this email directly, view it on GitHub [3], or unsubscribe [4]. You are receiving this because you are subscribed to this thread.Message ID: @.***>

--


Dr. Marinou Eleni @.*** IAASARS www.astro.noa.gr [5] National Observatory of Athens (NOA) I. Metaxa & V. Pavlou, P. Penteli, 15236, Athens

Remote sensing of Aerosols, Clouds & Trace gases Group (ReACT; https://react.space.noa.gr/) .... o o o o .. -\<, `\<, `\<, `\<, ...()/() ()/ () ()/ () ()/ ()


Links:

[1] https://github.com/thgeorgiou [2] https://github.com/HolgerPollyNet [3] https://github.com/NOA-ReACT/PollyXT-SCC-Pipelines/issues/26 [4] https://github.com/notifications/unsubscribe-auth/AQBHLQIHHFNOWMGRSVY7MYLWWHX5FANCNFSM6AAAAAAUTSSWO4 [5] http://www.astro.noa.gr

ulysses78 commented 1 year ago

As I mentioned earlier, the big drawback with this updated scripts of mine:

one has to add the 1064nm depol-channel to all sites in the locations.ini. In other words... if all the users of this repo updating their scripts, they also have to change their locations.ini file (adding pseudo 1064nm depol-channels). Otherwise the script will not work.

Of course I can make a new branch wih my updates, but I think it is best, if you change the code snippets in a way, that it is not necessary to change all the locations.ini entries (in terms of adding pseudo 1064nm depol-channels).

Do you agree?

thgeorgiou commented 1 year ago

Hello @ulysses78 ,

that is probably optimal. We will add your changes but make 1064nm optional. Do all Pollys (Pollies?) have both 355nm and 532nm channels? Do you think it is worth it to make all depol. channels optional?

Cheers

ulysses78 commented 1 year ago

Holger mentioned a few days ago, that not all pollys have two depol-channels. So yes, it would be best if all depol channels are optional. @HolgerPollyNet ... please correct me if I am wrong!

HolgerPollyNet commented 1 year ago

@ulysses78 Correct!

ulysses78 commented 1 year ago

Hi @elmarinou , hi @EvgeniaPapak, is there any progress in this issue so far?

Best, Andi

thgeorgiou commented 1 year ago

Hello @ulysses78, I have pushed version 1.13.0 that should cover this use-case. Please do let us know if everything works correctly!

ulysses78 commented 1 year ago

Hello @thgeorgiou , thank you! I upgraded pollyxt-pipelines via pip to 1.13.0. At the moment the script doesn't recognize the 1064nm channel, which I have added to the locations.ini file. It only generates the 355 and 532 nm vol-depol-channels. I found the new file pollyxt_pipelines/locations/__init__1064depol.py. Do I have to implement this one somehow, somewhere?

ulysses78 commented 1 year ago

@EvgeniaPapak .... do you have any idea? Am I doing something wrong with that? Do I have to flip a switch somewhere to trigger that 1064nm depolchannel-reading?

thgeorgiou commented 1 year ago

@ulysses78 Hello! I apologize for the late response, I missed the notification. Could I have a sample file with all channels to use for debugging?

ulysses78 commented 1 year ago

You can download some files from here:

https://owncloud.gwdg.de/index.php/f/2322020857

Can you access it?

Fuhrtermore I have attached the entry in our locations.ini file:

[CVO]
scc_code = cvo
lat = 16.87
lon = -24.99
altitude_asl = 10
daytime_configuration = 712
nighttime_configuration = 713
channel_id = 1995,1996,1997,1998,1999,2000,2001,2002,2005,2006,2007,2008,2009,2004,2003
background_low = 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0
background_high = 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240, 240
lr_input = 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1
temperature = 15
pressure = 1000
total_channel_355_nm_idx = 0
cross_channel_355_nm_idx = 1
total_channel_532_nm_idx = 4
cross_channel_532_nm_idx = 5
total_channel_1064_nm_idx = 7
cross_channel_1064_nm_idx = 14
calibration_355nm_total_channel_ids = 2010, 2011
calibration_355nm_cross_channel_ids = 2012, 2013
calibration_532nm_total_channel_ids = 2014, 2015
calibration_532nm_cross_channel_ids = 2016, 2017
calibration_1064nm_total_channel_ids = 2018, 2019
calibration_1064nm_cross_channel_ids = 2020, 2021
calibration_configuration_355nm = 714
calibration_configuration_532nm = 715
calibration_configuration_1064nm = 716
depol_calibration_zero_state = 0
profile_name = mindelo
sounding_provider = cloudnet
ulysses78 commented 1 year ago

Some more Info: This is the only entry/station in our locations.ini file with a 1064nm depolchannel. All the other stations do not have a 1064nm channel.

thgeorgiou commented 1 year ago

@ulysses78 Unfortunately I cannot access that link, I get a login prompt. Thanks for the location.ini, would it be OK with you if we added it to the release?

ulysses78 commented 1 year ago

This link should work now:

https://owncloud.gwdg.de/index.php/s/qTzwAFNDeHThmhV

ulysses78 commented 1 year ago

@ulysses78 Unfortunately I cannot access that link, I get a login prompt. Thanks for the location.ini, would it be OK with you if we added it to the release?

Yes, you can add it to the release.

thgeorgiou commented 1 year ago

Thank you for the files, I will come back to you later today with any fixes.

thgeorgiou commented 1 year ago

image

@ulysses78 Hello again! I can't reproduce the problem, the program creates 1064 calibration files as expected. Could you re-install the package just to be safe?

pip install  --force-reinstall pollyxt_pipelines

Or if you use pipx:

pipx reinstall pollyxt_pipelines

Thanks!

ulysses78 commented 1 year ago

I reinstalled pollyxt_pipelines, but it is not working for me. I even switched to another user, ceated a new python environment and installed pollyxt_pipelines as brand new. Same problem... only 355 and 532nm depolchannel files are created. And I am sure, that the correct locations.ini file is used, because I had to copy/paste it in the new environments locations.ini file. Otherwise the CVO location wouldn't be available.

thgeorgiou commented 1 year ago

@ulysses78 Thank you for looking into it, I think the pypi release was broken for 1.13.0. Could you please try again w/ the new release? (1.13.2)

ulysses78 commented 1 year ago

Yes, it is working now. Thanks a lot!

thgeorgiou commented 1 year ago

Great! Cheers