credit.interp ============= .. py:module:: credit.interp .. autoapi-nested-parse:: Functions to interpolate data to pressure and height coordinates. Functions --------- .. autoapisummary:: credit.interp.full_state_pressure_interpolation credit.interp.fast_state_interp_loop credit.interp.create_pressure_grid credit.interp.create_reduced_pressure_grid credit.interp.geopotential_from_model_vars credit.interp.interp_hybrid_to_pressure_levels credit.interp.interp_pressure_to_hybrid_levels credit.interp.interp_hybrid_to_hybrid_levels credit.interp.interp_geopotential_to_pressure_levels credit.interp.interp_temperature_to_pressure_levels credit.interp.interp_hybrid_to_height_agl credit.interp.mean_sea_level_pressure credit.interp.mean_sea_level_pressure_simple Module Contents --------------- .. py:function:: full_state_pressure_interpolation(state_dataset: xarray.Dataset, surface_geopotential: numpy.ndarray, pressure_levels: numpy.ndarray = np.array([500.0, 850.0]), interp_fields: tuple[str] = ('U', 'V', 'T', 'Q'), pres_ending: str = '_PRES', height_levels: numpy.ndarray = None, height_ending: str = '_HEIGHT', temperature_var: str = 'T', q_var: str = 'Q', surface_pressure_var: str = 'SP', geopotential_var: str = 'Z', time_var: str = 'time', lat_var: str = 'latitude', lon_var: str = 'longitude', pres_var: str = 'pressure', level_var: str = 'level', height_var: str = 'height_agl', model_level_file: str = '../credit/metadata/ERA5_Lev_Info.nc', verbose: int = 1, a_model_name: str = 'a_model', b_model_name: str = 'b_model', a_half_name: str = 'a_half', b_half_name: str = 'b_half', P0: float = 1.0, pressure_3d_var: str = 'P', mslp_temp_height: float = 1000.0, use_simple_mslp: bool = False, surface_geopotential_var: str = 'geopotential_at_surface') -> xarray.Dataset Interpolate the full state of the model to pressure and height coordinates. Interpolate full model state variables from model levels to pressure levels and height levels. The raw CREDIT model output are on hybrid sigma-pressure vertical levels, which start as terrain following near the surface and relax to constant pressure levels aloft. The state variables for CREDIT models (and hydrostatic models more generally) are u, v, temperature, specific humidity, and surface pressure, with surface geopotential as a static variable. To perform pressure and height interpolation, the following steps happen: 1. Pressure is calculated on every full model level (middle of the vertical grid box) and half level (top and bottom of the vertical grid box starting at the surface and ending at the model top of the atmospshere). This requires knowing the a and b coefficients for the model levels. 2. Geopotential on each hybrid sigma-pressure level is calculated from surface geopotential, pressure, temperature, and specific humidity as a vertical integral calculation. The calculation is sensitive to numerical precision, so data are cast to float64 before calling the geopotential calculation. 3. Interpolation is done from hybrid sigma-pressure levels to fixed pressure levels. Everything is interpolated linearly with log(pressure) as the x coordinate. For pressure levels below ground (where pressure of level > surface pressure), special extrapolation routines are done for temperature and geopotential while constant extrapolation is assumed for u, v, and q. 4. Interpolation to height above ground level is also performed at the end. Heights are defined in meters. :param state_dataset: state variables being interpolated :type state_dataset: xr.Dataset :param surface_geopotential: surface geopotential levels in units m^2/s^2. :type surface_geopotential: np.ndarray :param pressure_levels: pressure levels for interpolation in hPa. :type pressure_levels: np.ndarray :param interp_fields: fields to be interpolated. :type interp_fields: tuple[str] :param pres_ending: ending string to attach to pressure interpolated variables. :type pres_ending: str :param height_levels: height levels for interpolation to height above ground level in meters. :type height_levels: np.ndarray :param height_ending: ending string to attach to height interpolated variables. :type height_ending: str :param temperature_var: temperature variable to be interpolated (units K). :type temperature_var: str :param q_var: mixing ratio/specific humidity variable to be interpolated (units kg/kg). :type q_var: str :param surface_pressure_var: surface pressure variable (units Pa). :type surface_pressure_var: str :param geopotential_var: geopotential variable being derived (units m^2/s^2). :type geopotential_var: str :param time_var: time coordinate :type time_var: str :param lat_var: latitude coordinate :type lat_var: str :param lon_var: longitude coordinate :type lon_var: str :param pres_var: pressure coordinate :type pres_var: str :param level_var: name of level coordinate :type level_var: str :param height_var: height coordinate :type height_var: str :param model_level_file: relative path to file containing model levels. :type model_level_file: str :param verbose: verbosity level. If verbose > 0, print progress. :type verbose: int :param a_model_name: Name of A weight at level midpoints in sigma coordinate formula. 'a_model' by default. :type a_model_name: str :param b_model_name: Name of B weight at level midpoints in sigma coordinate formula. 'b_model' by default. :type b_model_name: str :param a_half_name: Name of A weight at level interfaces in sigma coordinate formula. 'a_half' by default. :type a_half_name: str :param b_half_name: Name of B weight at level interfaces in sigma coordinate formula. 'b_half' by default. :type b_half_name: str :param P0: reference pressure if pressure needs to be scaled. :type P0: float :param pressure_3d_var: Name of the 3D pressure field derived on the model grid. :type pressure_3d_var: str :param mslp_temp_height: height above ground level in meters where temperature is sampled for mslp calculation. :type mslp_temp_height: float :param use_simple_mslp: Whether to use the simple or complex MSLP calculation. :type use_simple_mslp: bool :param surface_geopotential_var: Name of the unnormalized surface geopotential variable being used for MSLP calculations :type surface_geopotential_var: str :returns: Dataset containing pressure interpolated variables. :rtype: pressure_ds (xr.Dataset) .. py:function:: fast_state_interp_loop(surface_pressure_data, state_dict, surface_geopotential, temperature_var, q_var, interp_fields, geopotential_var, a_model, b_model, a_half_full, b_half_full, pressure_levels, pres_ending, height_ending, height_levels, pressure_3d_var, level_var, levels) .. py:function:: create_pressure_grid(surface_pressure, model_a_half, model_b_half) Create a pressure 3D grid from a full set of vertical levels. Create a 3D pressure field at model levels from the surface pressure field and the hybrid sigma-pressure coefficients from ECMWF. Conversion is `pressure_3d = a + b * SP`. :param surface_pressure: (time, latitude, longitude) or (latitude, longitude) grid in units of Pa. :type surface_pressure: np.ndarray :param model_a_half: a coefficients at each model level being used in units of Pa. :type model_a_half: np.ndarray :param model_b_half: b coefficients at each model level being used (unitness). :type model_b_half: np.ndarray :returns: 3D pressure field with dimensions of surface_pressure and number of levels from model_a and model_b. :rtype: pressure_3d .. py:function:: create_reduced_pressure_grid(surface_pressure, model_a_full, model_b_full) Create a pressure 3D grid using sparse vertical levels. Create a 3D pressure field at model levels from the surface pressure field and the reduced set of hybrid sigma- pressure levels used in the CREDIT models. This function assumes that the coefficients for the full levels are being passed and then derives the half levels by taking the geometric means of the a and b coefficients on full levels. Conversion is `pressure_3d = a + b * SP`. :param surface_pressure: (time, latitude, longitude) or (latitude, longitude) grid in units of Pa. :type surface_pressure: np.ndarray :param model_a_full: a coefficients at each model level being used in units of Pa. :type model_a_full: np.ndarray :param model_b_full: b coefficients at each model level being used (unitless). :type model_b_full: np.ndarray :returns: 3D pressure field with dimensions of surface_pressure and number of levels from model_a and model_b. :rtype: pressure_3d .. py:function:: geopotential_from_model_vars(surface_geopotential, surface_pressure, temperature, specific_humidity, half_pressure) Calculate geopotential from model level data. Calculate geopotential from the base state variables. Geopotential height is calculated by adding thicknesses calculated within each half-model-level to account for variations in temperature and moisture between grid cells. Note that this function is calculating geopotential in units of (m^2 s^-2) not geopential height. To convert geopotential to geopotential height, divide geopotential by g (9.806 m s^-2). Geopotential height is defined as the height above mean sea level. To get height above ground level, substract the surface geoptential height field from the 3D geopotential height field. :param surface_geopotential: Surface geopotential in shape (y,x) and units m^2 s^-2. :type surface_geopotential: np.ndarray :param surface_pressure: Surface pressure in shape (y, x) and units Pa :type surface_pressure: np.ndarray :param temperature: temperature in shape (levels, y, x) and units K :type temperature: np.ndarray :param specific_humidity: mixing ratio in shape (levels, y, x) and units kg/kg. :type specific_humidity: np.ndarray :returns: geopotential on model levels in shape (levels, y, x) :rtype: model_geoptential (np.ndarray) .. py:function:: interp_hybrid_to_pressure_levels(model_var, model_pressure, interp_pressures, use_log=True) Interpolate to pressure levels. Interpolate data field from hybrid sigma-pressure vertical coordinates to pressure levels. `model_pressure` and `interp_pressure` should have consistent units with each other. :param model_var: 3D field on hybrid sigma-pressure levels with shape (levels, y, x). :type model_var: np.ndarray :param model_pressure: 3D pressure field with shape (levels, y, x) in units Pa or hPa :type model_pressure: np.ndarray :param interp_pressures: (np.ndarray): pressure levels for interpolation in units Pa or hPa. :param use_log: If True, use the natural logarithm of the pressure as the interpolation coordinate. Otherwise, use the pressure. :type use_log: bool :returns: 3D field on pressure levels with shape (len(interp_pressures), y, x). :rtype: pressure_var (np.ndarray) .. py:function:: interp_pressure_to_hybrid_levels(pressure_var, pressure_levels, model_pressure, surface_pressure) Interpolate fields on pressure levels to hybrid levels. Interpolate data field from hybrid sigma-pressure vertical coordinates to pressure levels. `model_pressure` and `pressure_levels` and 'surface_pressure' should have consistent units with each other. :param pressure_var: 3D field on pressure levels with shape (levels, y, x). :type pressure_var: np.ndarray :param pressure_levels: pressure levels for interpolation in units Pa or hPa. :type pressure_levels: np.double :param model_pressure: 3D pressure field with shape (levels, y, x) in units Pa or hPa :type model_pressure: np.ndarray :param surface_pressure: pressure at the surface in units Pa or hPa. :type surface_pressure: np.ndarray :returns: 3D field on hybrid sigma-pressure levels with shape (model_pressure.shape[0], y, x). :rtype: model_var (np.ndarray) .. py:function:: interp_hybrid_to_hybrid_levels(hybrid_var, hybrid_pressure, target_pressure) Interpolate fields on hybrid levels to hybrid levels via pressure. Interpolate data from hybrid sigma-pressure vertical coordinates to other hybrid levels. :param hybrid_var: 3D field on hybrid sigma-pressure levels with shape (levels, y, x). :type hybrid_var: np.ndarray :param hybrid_pressure: pressure levels for interpolation in units Pa or hPa. :type hybrid_pressure: np.double :param target_pressure: 3D target pressure fields with shape (levels, y, x) in units Pa or hPa :type target_pressure: np.ndarray :returns: 3D field on hybrid sigma-pressure levels with shape (target_pressure.shape[0], y, x). :rtype: model_var (np.ndarray) .. py:function:: interp_geopotential_to_pressure_levels(geopotential, model_pressure, interp_pressures, surface_pressure, surface_geopotential, temperature_k, temp_height=150) Interpolate geopotential field to pressure levels. Interpolate geopotential field from hybrid sigma-pressure vertical coordinates to pressure levels. `model_pressure` and `interp_pressure` should have consistent units of hPa or Pa. Geopotential height is extrapolated below the surface based on Eq. 15 in Trenberth et al. (1993). :param geopotential: geopotential in units m^2/s^2. :type geopotential: np.ndarray :param model_pressure: 3D pressure field with shape (levels, y, x) in units Pa or hPa :type model_pressure: np.ndarray :param interp_pressures: pressure levels for interpolation in units Pa or hPa. :type interp_pressures: np.ndarray :param surface_pressure: pressure at the surface in units Pa or hPa. :type surface_pressure: np.ndarray :param surface_geopotential: geopotential at the surface in units m^2/s^2. :type surface_geopotential: np.ndarray :param temperature_k: temperature in units K. :type temperature_k: np.ndarray :param temp_height: height above ground of nearest vertical grid cell. :type temp_height: float :returns: 3D field on pressure levels with shape (len(interp_pressures), y, x). :rtype: pressure_var (np.ndarray) .. py:function:: interp_temperature_to_pressure_levels(model_var, model_pressure, interp_pressures, surface_pressure, surface_geopotential, geopotential, temp_height=150) Interpolate the temperature field to pressure levels. Interpolate temperature field from hybrid sigma-pressure vertical coordinates to pressure levels. `model_pressure` and `interp_pressure` should have consistent units of hPa or Pa. Temperature is extrapolated below the surface based on Eq. 16 in Trenberth et al. (1993). :param model_var: 3D field on hybrid sigma-pressure levels with shape (levels, y, x). :type model_var: np.ndarray :param model_pressure: 3D pressure field with shape (levels, y, x) in units Pa :type model_pressure: np.ndarray :param interp_pressures: (np.ndarray): pressure levels for interpolation in units Pa or. :param surface_pressure: pressure at the surface in units Pa or hPa. :type surface_pressure: np.ndarray :param surface_geopotential: geopotential at the surface in units m^2/s^2. :type surface_geopotential: np.ndarray :param temp_height: height above ground of nearest vertical grid cell. :type temp_height: float :returns: 3D field on pressure levels with shape (len(interp_pressures), y, x). :rtype: pressure_var (np.ndarray) .. py:function:: interp_hybrid_to_height_agl(model_var: numpy.ndarray, interp_heights_m: numpy.ndarray, geopotential: numpy.ndarray, surface_geopotential: numpy.ndarray) Interpolate data on hybrid sigma-pressure levels to heights above ground level in meters. :param model_var: State variable of shape [levels, lat, lon] :type model_var: np.ndarray :param interp_heights_m: 1D array of height levels in meters above ground level. :type interp_heights_m: np.ndarray :param geopotential: geopotential on model levels in units of m^2/s^2. :type geopotential: np.ndarray :param surface_geopotential: geopotential at the surface in units of m^2/s^2. :returns: State variable on height above ground levels in shape [interp_heights, lat, lon]. :rtype: height_var (np.ndarray) .. py:function:: mean_sea_level_pressure(surface_pressure_pa, temperature_k, pressure_pa, surface_geopotential, geopotential, temp_height=150.0) Calculate the mean sea level pressure. Calculate mean sea level pressure from surface pressure, lowest model level temperature, the pressure of the lowest model level (derived from create_pressure_grid), and surface_geopotential. This calculation is based on the procedure from Trenberth et al. (1993) implemented in CESM CAM. Trenberth, K., J. Berry , and L. Buja, 1993: Vertical Interpolation and Truncation of Model-Coordinate, University Corporation for Atmospheric Research, https://doi.org/10.5065/D6HX19NH. CAM implementation: https://github.com/ESCOMP/CAM/blob/cam_cesm2_2_rel/src/physics/cam/cpslec.F90 :param surface_pressure_pa: surface pressure in Pascals :param temperature_k: Temperature at the lowest model level in Kelvin. :param pressure_pa: Pressure at the lowest model level in Pascals. :param surface_geopotential: Geopotential of the surface in m^2 s^-2. :param geopotential: Geopotential at all levels. :param temp_height: height of nearest vertical grid cell :returns: Mean sea level pressure in Pascals. :rtype: mslp .. py:function:: mean_sea_level_pressure_simple(surface_pressure_pa, temperature_k, surface_geopotential) Simpler calculation for mean sea level pressure that only requires 2D fields of pressure (Pa), temperature (K), and surface geopotential (m ** 2 s ** -2). Based on Trenberth et al. 1993 calculation but simplified by removing the T* calculation since it seemed to only vary by about 0.2 K and requires a lot more data to compute. Trenberth, K., J. Berry , and L. Buja, 1993: Vertical Interpolation and Truncation of Model-Coordinate, University Corporation for Atmospheric Research, https://doi.org/10.5065/D6HX19NH. :param surface_pressure_pa: surface pressure in Pascals :param temperature_k: temperature in Kelvin :param surface_geopotential: surface geopotential in m^2 s^-2. If you have surface height, multiply by g (9.81 m2s-2) :returns: mean sea level pressure in Pascals.