emc2.core (emc2.core)¶
The procedures in this module contain the core data structures of EMC^2. In particular, the Instrument class that describes the characteristics of the instrument to be simulated is stored here. In addition, global constants used by EMC^2 are also stored in this module.
|
This stores the information for the KAZR (Ka-band radar). |
|
This stores the information for 532 nm lidars ,e.g., the High Spectral Resolution Lidar (HSRL), micropulse lidar (MPL). |
|
This stores the information for the ARM CSAPR. |
|
This stores the information for the NOAA NEXRAD radar Based on https://www.roc.noaa.gov/WSR88D/Engineering/NEXRADTechInfo.aspx. |
|
This stores the information for the AMF XSACR. |
|
This stores the information for the 1064 nm lidars, e.g., the 2-ch HSRL. |
|
This loads a ModelE simulation with all of the necessary parameters for EMC^2 to run. |
This is a test Model structure used only for unit testing. |
|
|
This loads an E3SM simulation output with all of the necessary parameters for EMC^2 to run. |
|
This loads a DHARMA simulation with all of the necessary parameters for EMC^2 to run. |
|
This load a WRF simulation and all of the necessary parameters from the simulation. |
This is a test Model structure used only for unit testing. |
|
This is a test Model structure used only for unit testing. |
|
This is a test Model structure used only for unit testing. |
|
This is a test Model structure used only for unit testing. |
|
|
This is the base class which holds the information needed to contain the instrument parameters for the simulator. |
|
This class stores the model specific parameters for the radar simulator. |
In addition the emc2.core.quantity()
is equivalent to pint’s
UnitRegistry().Quantity object. This allows for the use of quantities
with units, making EMC^2 unit aware.
- class emc2.core.Model[source]¶
This class stores the model specific parameters for the radar simulator.
- Attributes:
- Rho_hyd: dict
A dictionary whose keys are the names of the model’s hydrometeor classes and whose values are the density of said hydrometeors in \(kg\ m^{-3}\)
- fluffy: dict
A dictionary whose keys are the names of the model’s ice hydrometeor classes and whose values are the ice fluffiness factor for the fwd calculations using r_e, where values of 0 - equal volume sphere, 1 - fluffy sphere i.e., diameter = maximum dimension.
- lidar_ratio: dict
A dictionary whose keys are the names of the model’s hydrometeor classes and whose values are the lidar_ratio of said hydrometeors.
- vel_param_a: dict
A dictionary whose keys are the names of the model’s hydrometeor classes and whose values are the \(a\) parameters to the equation \(V = aD^b\) used to calculate terminal velocity corresponding to each hydrometeor.
- vel_param_b: dict
A dictionary whose keys are the names of the model’s hydrometeor classes and whose values are the \(b\) parameters to the equation \(V = aD^b\) used to calculate terminal velocity corresponding to each hydrometeor.
- N_field: dict
A dictionary whose keys are the names of the model’s hydrometeor classes and whose values are the number concentrations in \(cm^{-3}\) corresponding to each hydrometeor class.
- T_field: str
A string containing the name of the temperature field in the model.
- q_field: str
A string containing the name of the water vapor mixing ratio field (in kg/kg) in the model.
- p_field: str
A string containing the name of the pressure field (in mbar) in the model.
- z_field: str
A string containing the name of the height field (in m) in the model.
- conv_frac_names: dict
A dictionary containing the names of the convective fraction corresponding to each hydrometeor class in the model.
- strat_frac_names: dict
A dictionary containing the names of the stratiform fraction corresponding to each hydrometeor class in the model.
- conv_frac_names_for_rad: dict
A dictionary containing the names of the convective fraction corresponding to each hydrometeor class in the model for the radiation scheme.
- strat_frac_names_for_rad: dict
A dictionary containing the names of the stratiform fraction corresponding to each hydrometeor class in the model for the radiation scheme.
- conv_re_fields: dict
A dictionary containing the names of the effective radii of each convective hydrometeor class
- strat_re_fields: dict
A dictionary containing the names of the effective radii of each stratiform hydrometeor class
- asp_ratio_func: dict
A dictionary that returns hydrometeor aspect ratios as a function of maximum dimension in mm.
- hyd_types: list
list of hydrometeor classes to include in calcuations. By default set to be consistent with the model represented by the Model subclass.
- time_dim: str
The name of the time dimension in the model.
- height_dim: str
The name of the height dimension in the model.
- lat_dim: str
Name of the latitude dimension in the model (relevant for regional output)
- lon_dim: str
Name of the longitude dimension in the model (relevant for regional output)
- mcphys_scheme: str
Name of the microphysics scheme to use for models with multiple microphysics schemes.
- stacked_time_dim: str or None
This attribute becomes a string of the original time dimension name only if stacking was required to enable EMC2 to processes a domain output (time x lat x lon).
- process_conv: bool
If True, then processing convective model output (can typically be set to False for some models).
- model_name: str
The name of the model.
- variable_density: dict
If the model allows for particle density for vary (e.g. 2-moment NSSL), then this is a dict pointing to the variable with the density for each hydrometeor class
- check_and_stack_time_lat_lon(out_coord_name='time_lat_lon', file_path=None, order_dim=True)[source]¶
Stack the time dim together with the lat and lon dims (if the lat and/or lon dims are longer than 1) to enable EMC^2 processing of regional model output. Otherwise, squeezing the lat and lon dims (if they exist in dataset). Finally, the method reorder dimensions to time x height for proper processing by calling the “permute_dims_for_processing” class method. NOTE: tmp variables for lat, lon, and time are produced as xr.Datasets still have many unresolved bugs associated with pandas multi-indexing implemented in xarray for stacking (e.g., new GitHub issue #5692). Practically, after the subcolumn processing the stacking information is lost so an alternative dedicated method is used for unstacking
- Parameters:
- out_coord_name: str
Name of output stacked coordinate.
- file_path: str
Path and filename of model simulation output.
- order_dim: bool
When True, reorder dimensions to time x height for proper processing.
- finalize_subcol_fields(more_fieldnames=[])[source]¶
Remove all zero values from subcolumn output fields enabling better visualization. Can be applied over additional fields using the more_fieldnames input parameter.
- property hydrometeor_classes¶
The list of hydrometeor classes.
- load_subcolumns_from_netcdf(file_name)[source]¶
Load all of the subcolumn data from a previously saved netCDF file. The dataset being loaded must match the current number of subcolumns if there are any generated.
- Parameters:
- file_name: str
Name of the file to save.
- property num_hydrometeor_classes¶
The number of hydrometeor classes used
- property num_subcolumns¶
Gets the number of subcolumns in the model. Will return 0 if the number of subcolumns has not yet been set.
- permute_dims_for_processing(base_order=None, base_dim_first=True)[source]¶
Reorder dims for consistent processing such that the order is: subcolumn x time x height. Note: lat/lon dims are assumed to already be stacked with the time dim.
- Parameters:
- base_order: list or None
List of preffered dimension order. Use default if None
- base_dim_first: bool
Make the base dims (height and time) the first ones in the permutation if True.
- remove_appended_str(all_appended_in_lat=False)[source]¶
Remove appended strings from xr.Dataset coords and fieldnames based on lat/lon coord names (typically required when using post-processed output data files).
- Parameters:
- all_appended_in_lat: bool
If True using only the appended str portion to the lat_dim. Otherwise, combining the appended str from both the lat and lon dims (relevant if appended_str is True).
- remove_subcol_fields(cloud_class='conv')[source]¶
Remove all subcolumn output fields for the given cloud class to save memory (mainly releveant for CESM and E3SM).
- subcolumns_to_netcdf(file_name)[source]¶
Saves all of the simulated subcolumn parameters to a netCDF file.
- Parameters:
- file_name: str
The name of the file to save to.
- unstack_time_lat_lon(order_dim=True, squeeze_single_dims=True)[source]¶
Unstack the time, lat, and lon dims if they were previously stacked together (self.stacked_time_dim is not None). Finally, the method reorder dimensions to time x height for proper processing by calling the “permute_dims_for_processing” class method. Note (relevant for squeeze_single_dims == True): If the time dimension size is 1, that dimension is squeezed. Similarly, if the number of subcolumns is 1 (subcolumn generator is turned off), the subcolumn dimension is squeezed as well. NOTE: This is a dedicated method written because xr.Datasets still have many unresolved bugs associated with pandas multi-indexing implemented in xarray for stacking (e.g., new GitHub issue #5692). Practically, after the subcolumn processing the stacking information is lost so this is an alternative dedicated method.
- Parameters:
- order_dim: bool
When True, reorder dimensions to subcolumn x time x height for proper processing.
- squeeze_single_dims: bool
If True, squeezing the time and/or subcolumn dimension(s) if their lengh is 1.
- class emc2.core.Instrument(frequency=None, wavelength=None)[source]¶
This is the base class which holds the information needed to contain the instrument parameters for the simulator.
- Attributes:
- instrument_str: str
The name of the instrument.
- instrument_class: str
The class of the instrument. Currently must be one of ‘radar,’ or ‘lidar’.
- freq: float
The frequency of the instrument.
- wavelength: float
The wavelength of the instrument
- beta_p_phase_thresh: list of dicts or None
If a list, each index contains a dictionaly with class name, class integer value (mask order), LDR value bounds, and the corresponding beta_p threshold (thresholds are linearly interpolated between LDR values). In order for the method to operate properly, the list should be arranged from the lowest to highest beta_p threshold values for a given LDR, that is, beta_p[i+1 | LDR=x] >= beta_p[i | LDR=x].
- ext_OD: float
The optical depth where we have full extinction of the lidar signal.
- OD_from_sfc: Bool
If True (default), optical depth will be calculated from the surface. If False, optical depth will be calculated from the top of the atmosphere.
- eta: float
Multiple scattering coefficient.
- K_w: float
The index of refraction of water used for Ze calculation. See the ARM KAZR handbook (Widener et al. 2012)
- eps_liq: float
The complex dielectric constant for liquid water.
- pt: float
Transmitting power in Watts.
- theta: float
3 dB beam width in degrees
- gain: float
The antenna gain in linear units.
- Z_min_1km: float
The minimum detectable signal at 1 km in dBZ
- lr: float
Attenuation based on the the general attributes in the spectra files.
- pr_noise_ge: float
Minimum detectable signal in mW.
- tau_ge: float
Pulse width in mus.
- tau_md: float
Pulse width in mus.