pycontrails.models.cocip.Cocip¶
- class pycontrails.models.cocip.Cocip(met, rad, params=None, **params_kwargs)¶
Bases:
ModelContrail Cirrus Prediction Model (CoCiP).
Published by Ulrich Schumann et. al. (DLR Institute of Atmospheric Physics) in [Schumann, 2012], [Schumann et al., 2012].
- Parameters:
met (
MetDataset) – Pressure level dataset containingmet_variablesvariables. See Notes for variable names by data source.rad (
MetDataset) – Single level dataset containing top of atmosphere radiation fluxes. See Notes for variable names by data source.params (
dict[str,Any] | None, optional) – Override Cocip model parameters with dictionary. SeeCocipFlightParamsfor model parameters.**params_kwargs (
Any) – Override Cocip model parameters with keyword arguments. SeeCocipFlightParamsfor model parameters.
Notes
Inputs
The required meteorology variables depend on the data source.
Cocipsupports data-source-specific variables from ECMWF models (HRES, ERA5) and the NCEP GFS, plus a generic set of model-agnostic variables.See
met_variablesandrad_variablesfor the list of required variables to themetandradparameters, respectively. When an item in one of these arrays is atuple, variable keys depend on data source.A warning will be raised if meteorology data is from a source not currently supported by a pycontrails datalib. In this case it is the responsibility of the user to ensure that meteorology data is formatted correctly. The warning can be suppressed with a context manager:
import warnings with warnings.catch_warnings(): warnings.simplefilter("ignore", category=UserWarning, message="Unknown provider") cocip = Cocip(met, rad, ...)
The current list of required variables (labelled by
"standard_name"):Variable keys for pressure level data¶ Parameter
ECMWF
GFS
Generic
Air Temperature
air_temperatureair_temperatureair_temperatureSpecific Humidity
specific_humidityspecific_humidityspecific_humidityEastward wind
eastward_windeastward_windeastward_windNorthward wind
northward_windnorthward_windnorthward_windVertical velocity
lagrangian_tendency_of_air_pressurelagrangian_tendency_of_air_pressurelagrangian_tendency_of_air_pressureIce water content
specific_cloud_ice_water_contentice_water_mixing_ratiomass_fraction_of_cloud_ice_in_airVariable keys for single-level radiation data¶ Parameter
ECMWF
GFS
Generic
Top solar radiation
top_net_solar_radiationtoa_upward_shortwave_fluxtoa_net_downward_shortwave_fluxTop thermal radiation
top_net_thermal_radiationtoa_upward_longwave_fluxtoa_outgoing_longwave_fluxModifications
This implementation differs from original CoCiP (Fortran) implementation in a few places:
This model uses aircraft performance and emissions models to calculate nvPM, fuel flow, and overall propulsion efficiency, if not already provided.
As described in [Teoh et al., 2022], this implementation sets the initial ice particle activation rate to be a function of the difference between the ambient temperature and the critical SAC threshold temperature. See
pycontrails.models.sac.T_critical_sac().Isobaric heat capacity calculation. The original model uses a constant value of 1004 \(J \ kg^{-1} \ K^{-1}\), whereas this model calculates isobaric heat capacity as a function of specific humidity. See
pycontrails.physics.thermo.c_pm().Solar direct radiation. The original algorithm uses ECMWF radiation variable tisr (top incident solar radiation) as solar direct radiation value. This implementation calculates the theoretical solar direct radiation at any arbitrary point in the atmosphere. See
pycontrails.physics.geo.solar_direct_radiation().Segment angle. The segment angle calculations for flights and contrail segments have been updated to use more precise spherical geometry instead of a triangular approximation. As the triangle approaches zero, the two calculations agree. See
pycontrails.physics.geo.segment_angle().Integration. This implementation consistently uses left-Riemann sums in the time integration of contrail segments.
Segment length ratio. Instead of taking a geometric mean between contrail segments before/after advection, a simple ratio is computed. See
contrail_properties.segment_length_ratio().Segment energy flux. This implementation does not average spatially contiguous contrail segments when calculating the mean energy flux for the segment of interest. See
contrail_properties.mean_energy_flux_per_m().
This implementation is regression tested against results from [Teoh et al., 2022].
Outputs
NaN values may appear in model output. Specifically,
np.nanvalues are used to indicate:Flight waypoint or contrail waypoint is not contained with the
metdomain.The variable was NOT computed during the model evaluation. For example, at flight waypoints not producing any persistent contrails, “radiative” variables (
rsr,olr,rf_sw,rf_lw,rf_net) are not computed. Consequently, the corresponding values in the output ofeval()are NaN. One exception to this rule is found onef(energy forcing) contrail_age predictions. For these two “cumulative” variables, waypoints not producing any persistent contrails are assigned 0 values.
References
See also
CocipFlightParams,wake_vortex,contrail_properties,radiative_forcing,humidity_scaling,Emissions,sac,tau_cirrus- __init__(met, rad, params=None, **params_kwargs)¶
Methods
__init__(met, rad[, params])Return an ECMWF-specific list of required meteorology variables.
Return an ECMWF-specific list of required radiation variables.
eval([source])Run CoCiP simulation on flight.
Return a model-agnostic list of required meteorology variables.
Return a model-agnostic list of required radiation variables.
get_data_param(other, key[, default, set_attr])Get data from other source-compatible object with default set by model parameter key.
get_source_param(key[, default, set_attr])Get source data with default set by parameter key.
Return a GFS-specific list of required meteorology variables.
Return a GFS-specific list of required radiation variables.
Ensure that
metis a MetDataset.require_source_type(type_)Ensure that
sourceistype_.set_source([source])Attach original or copy of input
sourcetosource.set_source_met([optional, variable])Ensure or interpolate each required
met_variablesonsource.transfer_met_source_attrs([source])Transfer met source metadata from
mettosource.update_params([params])Update model parameters on
params.Attributes
Contrail evolution output from model.
xr.Datasetrepresentation of contrail evolution.List of
GeoVectorDatasetcontrail objects - one for each timestepRadiation data formatted as a
MetDatasetat a single pressure level [-1]Array of
numpy.datetime64time steps for contrail evolutionGenerate a unique hash for model instance.
Shortcut to create interpolation arguments from
params.Met data is not optional
Require meteorology is not None on __init__()
Required meteorology pressure level variables.
Additional met variables used to support outputs
Instantiated model parameters, in dictionary form
Minimal set of met variables needed to run the model after pre-processing.
Required single-level top of atmosphere radiation variables.
Last Flight modeled in
eval()- contrail¶
Contrail evolution output from model.
Set to None when no contrails are formed. Otherwise, this is a
pandas.DataFramedescribing the evolution of the contrail. Columns include:waypoint: The index of the waypoint in the original flight creating the contrail. This can be used to join the contrail DataFrame to thesource.formation_time: Time of contrail formation. Agrees with thetimecolumn insource.continuous: Boolean indicating whether the contrail is continuous or not.persistent: Boolean indicating whether the contrail is persistent or not. A contrail segment is considered continuous if both the current and the next contrail waypoint at the same time step persist.segment_length: Length of the contrail segment, [\(m\)].sin_a,cos_a: Sine and cosine of the segment angle.width,depth: Contrail width and depth, [\(m\)].sigma_yz: Theyzcomponent of the covariance matrix, [\(m^{2}\)]. Seecontrail_properties.plume_temporal_evolution().q_sat: Saturation specific humidity over ice, [\(kg \ kg^{-1}\)].n_ice_per_m: Number of ice particles per distance, [\(m^{-1}\)].iwc: Ice water content, [\(kg_{ice} kg_{air}^{-1}\)].tau_contrail: Optical depth of the contrail. Seecontrail_properties.contrail_optical_depth().rf_sw,rf_lw,rf_net: Shortwave, longwave, and net instantaneous radiative forcing, [\(W \ m^{-2}\)] at the contrail waypoint.ef: Energy forcing, [\(J\)] at the contrail waypoint. Seecontrail_properties.energy_forcing().
- contrail_dataset¶
xr.Datasetrepresentation of contrail evolution.
- contrail_list¶
List of
GeoVectorDatasetcontrail objects - one for each timestep
- default_params¶
alias of
CocipFlightParams
- downselect_met()¶
Downselect
metdomain to the max/min bounds ofsource.Override this method if special handling is needed in met down-selection.
sourcemust be defined before callingdownselect_met().This method copies and re-assigns
metusingmet.copy()to avoid side-effects.
- Raises:
ValueError – Raised if
sourceis not defined. Raised ifsourceis not aGeoVectorDataset.
- classmethod ecmwf_met_variables()¶
Return an ECMWF-specific list of required meteorology variables.
- Returns:
tuple[MetVariable]– List of ECMWF-specific variants of required variables
- classmethod ecmwf_rad_variables()¶
Return an ECMWF-specific list of required radiation variables.
- Returns:
tuple[MetVariable,]– List of ECMWF-specific variants of required variables
- eval(source=None, **params)¶
Run CoCiP simulation on flight.
Simulates the formation and evolution of contrails from a Flight using the contrail cirrus prediction model (CoCiP) from Schumann (2012) [Schumann, 2012].
Changed in version 0.25.11: Previously, any waypoint not surviving the wake vortex downwash phase of CoCiP was assigned a nan-value in the
efarray within the model output. This is no longer the case. Instead, energy forcing is set to 0.0 for all waypoints which fail to produce persistent contrails. In particular, nan values in theefarray are only used to indicate an out-of-met-domain waypoint. The same convention is now used for output variablescontrail_ageandcocipas well.- Parameters:
source (
Flight | Sequence[Flight] | None) – Input Flight(s) to model.**params (
Any) – Overwrite model parameters before eval.
- Returns:
Flight | list[Flight]– Flight(s) with updated Contrail data. The model parameter “verbose_outputs” determines the variables on the return flight object.
References
- classmethod generic_met_variables()¶
Return a model-agnostic list of required meteorology variables.
- Returns:
tuple[MetVariable]– List of model-agnostic variants of required variables
- classmethod generic_rad_variables()¶
Return a model-agnostic list of required radiation variables.
- Returns:
tuple[MetVariable,]– List of model-agnostic variants of required variables
- get_data_param(other, key, default=<object object>, *, set_attr=True)¶
Get data from other source-compatible object with default set by model parameter key.
Retrieves data with the following hierarchy:
other.data[key]. Returnsnp.ndarray | xr.DataArray.other.attrs[key]params[key]default
In case 3., the value of
params[key]is attached toother.attrs[key]unlessset_attris set to False.- Parameters:
- Returns:
Any– Value(s) found for key inotherdata,otherattrs, or model params- Raises:
KeyError – Raises KeyError if key is not found in any location and
defaultis not provided.
- get_source_param(key, default=<object object>, *, set_attr=True)¶
Get source data with default set by parameter key.
Retrieves data with the following hierarchy:
source.data[key]. Returnsnp.ndarray | xr.DataArray.source.attrs[key]params[key]default
In case 3., the value of
params[key]is attached tosource.attrs[key]unlessset_attris set to False.- Parameters:
- Returns:
Any– Value(s) found for key in source data, source attrs, or model params- Raises:
KeyError – Raises KeyError if key is not found in any location and
defaultis not provided.
- classmethod gfs_met_variables()¶
Return a GFS-specific list of required meteorology variables.
- Returns:
tuple[MetVariable]– List of GFS-specific variants of required variables
- classmethod gfs_rad_variables()¶
Return a GFS-specific list of required radiation variables.
- Returns:
tuple[MetVariable,]– List of GFS-specific variants of required variables
- property hash¶
Generate a unique hash for model instance.
- Returns:
str– Unique hash for model instance (sha1)
- property interp_kwargs¶
Shortcut to create interpolation arguments from
params.The output of this is useful for passing to
interpolate_met().- Returns:
dict[str,Any]– Dictionary with keys”method”
”bounds_error”
”fill_value”
”localize”
”use_indices”
”q_method”
as determined by
params.
- long_name = 'Contrail Cirrus Prediction Model'¶
- met¶
Met data is not optional
- met_required = True¶
Require meteorology is not None on __init__()
- met_variables = (MetVariable(short_name='t', standard_name='air_temperature', long_name='Air Temperature', level_type='isobaricInhPa', ecmwf_id=130, grib1_id=11, grib2_id=(0, 0, 0), units='K', amip='ta', description='Air temperature is the bulk temperature of the air, not the surface (skin) temperature.'), MetVariable(short_name='q', standard_name='specific_humidity', long_name='Specific Humidity', level_type='isobaricInhPa', ecmwf_id=133, grib1_id=51, grib2_id=(0, 1, 0), units='kg kg**-1', amip='hus', description='Specific means per unit mass. Specific humidity is the mass fraction of water vapor in (moist) air.'), MetVariable(short_name='u', standard_name='eastward_wind', long_name='Eastward Wind', level_type='isobaricInhPa', ecmwf_id=131, grib1_id=33, grib2_id=(0, 2, 2), units='m s**-1', amip='ua', description='"Eastward" indicates a vector component which is positive when directed eastward (negative westward). Wind is defined as a two-dimensional (horizontal) air velocity vector, with no vertical component.'), MetVariable(short_name='v', standard_name='northward_wind', long_name='Northward Wind', level_type='isobaricInhPa', ecmwf_id=132, grib1_id=34, grib2_id=(0, 2, 3), units='m s**-1', amip='va', description='"Northward" indicates a vector component which is positive when directed northward (negative southward). Wind is defined as a two-dimensional (horizontal) air velocity vector, with no vertical component.'), MetVariable(short_name='w', standard_name='lagrangian_tendency_of_air_pressure', long_name='Vertical Velocity (omega)', level_type='isobaricInhPa', ecmwf_id=135, grib1_id=39, grib2_id=(0, 2, 8), units='Pa s**-1', amip='wap', description='The Lagrangian tendency of air pressure, often called "omega", plays the role of the upward component of air velocity when air pressure is being used as the vertical coordinate. If the vertical air velocity is upwards, it is negative when expressed as a tendency of air pressure; downwards is positive. Air pressure is the force per unit area which would be exerted when the moving gas molecules of which the air is composed strike a theoretical surface of any orientation.'), (MetVariable(short_name='cli', standard_name='mass_fraction_of_cloud_ice_in_air', long_name='Mass fraction of cloud ice in air', level_type='isobaricInhPa', ecmwf_id=None, grib1_id=None, grib2_id=None, units='kg kg**-1', amip='cli', description='The mass fraction of cloud ice in moist air.'), MetVariable(short_name='ciwc', standard_name='specific_cloud_ice_water_content', long_name='Specific cloud ice water content', level_type='isobaricInhPa', ecmwf_id=247, grib1_id=None, grib2_id=(0, 1, 84), units='kg kg**-1', amip=None, description="This parameter is the mass of cloud ice particles per kilogram of the total mass of moist air. The 'total mass of moist air' is the sum of the dry air, water vapour, cloud liquid, cloud ice, rain and falling snow. This parameter represents the average value for a grid box."), MetVariable(short_name='icmr', standard_name='ice_water_mixing_ratio', long_name='Cloud ice water mixing ratio', level_type='isobaricInhPa', ecmwf_id=260019, grib1_id=None, grib2_id=(0, 1, 23), units='kg kg**-1', amip=None, description='This parameter is the mass of cloud ice particles per kilogram of the total mass of dry air. ')))¶
Required meteorology pressure level variables. Each element in the list is a
MetVariableor atuple[MetVariable]. If element is atuple[MetVariable], the variable depends on the data source and the tuple must include entries for a model-agnostic variable, an ECMWF-specific variable, and a GFS-specific variable. Only one of the three variable in the tuple is required for model evaluation.
- name = 'cocip'¶
- optional_met_variables = ((MetVariable(short_name='z', standard_name='geopotential', long_name='Geopotential', level_type='isobaricInhPa', ecmwf_id=129, grib1_id=6, grib2_id=(0, 3, 4), units='m**2 s**-2', amip=None, description='Geopotential is the sum of the specific gravitational potential energy relative to the geoid and the specific centripetal potential energy.'), MetVariable(short_name='gh', standard_name='geopotential_height', long_name='Geopotential Height', level_type='isobaricInhPa', ecmwf_id=156, grib1_id=7, grib2_id=(0, 3, 5), units='m', amip='zg', description='Geopotential is the sum of the specific gravitational potential energy relative to the geoid and the specific centripetal potential energy. Geopotential height is the geopotential divided by the standard acceleration due to gravity. It is numerically similar to the altitude (or geometric height) and not to the quantity with standard name height, which is relative to the surface.')), (MetVariable(short_name='cl', standard_name='cloud_area_fraction_in_atmosphere_layer', long_name='Cloud area fraction in atmosphere layer', level_type='isobaricInhPa', ecmwf_id=None, grib1_id=None, grib2_id=None, units='[0 - 1]', amip=None, description='The fraction of the horizontal area of a grid cell that contains cloud.'), MetVariable(short_name='cc', standard_name='fraction_of_cloud_cover', long_name='Cloud area fraction in atmosphere layer', level_type='isobaricInhPa', ecmwf_id=248, grib1_id=None, grib2_id=(0, 6, 32), units='[0 - 1]', amip='cl', description='This parameter is the proportion of a grid box covered by cloud (liquid or ice) at a specific pressure level.'), MetVariable(short_name='tcc', standard_name='total_cloud_cover_isobaric', long_name='Total cloud cover at isobaric surface', level_type='isobaricInhPa', ecmwf_id=228164, grib1_id=None, grib2_id=(0, 6, 1), units='%', amip=None, description='This parameter is the percentage of a grid box covered by cloud (liquid or ice) at a specific pressure level.')))¶
Additional met variables used to support outputs
Changed in version 0.48.0: Moved Geopotential from
met_variablestooptional_met_variables
- params¶
Instantiated model parameters, in dictionary form
- processed_met_variables = (MetVariable(short_name='t', standard_name='air_temperature', long_name='Air Temperature', level_type='isobaricInhPa', ecmwf_id=130, grib1_id=11, grib2_id=(0, 0, 0), units='K', amip='ta', description='Air temperature is the bulk temperature of the air, not the surface (skin) temperature.'), MetVariable(short_name='q', standard_name='specific_humidity', long_name='Specific Humidity', level_type='isobaricInhPa', ecmwf_id=133, grib1_id=51, grib2_id=(0, 1, 0), units='kg kg**-1', amip='hus', description='Specific means per unit mass. Specific humidity is the mass fraction of water vapor in (moist) air.'), MetVariable(short_name='u', standard_name='eastward_wind', long_name='Eastward Wind', level_type='isobaricInhPa', ecmwf_id=131, grib1_id=33, grib2_id=(0, 2, 2), units='m s**-1', amip='ua', description='"Eastward" indicates a vector component which is positive when directed eastward (negative westward). Wind is defined as a two-dimensional (horizontal) air velocity vector, with no vertical component.'), MetVariable(short_name='v', standard_name='northward_wind', long_name='Northward Wind', level_type='isobaricInhPa', ecmwf_id=132, grib1_id=34, grib2_id=(0, 2, 3), units='m s**-1', amip='va', description='"Northward" indicates a vector component which is positive when directed northward (negative southward). Wind is defined as a two-dimensional (horizontal) air velocity vector, with no vertical component.'), MetVariable(short_name='w', standard_name='lagrangian_tendency_of_air_pressure', long_name='Vertical Velocity (omega)', level_type='isobaricInhPa', ecmwf_id=135, grib1_id=39, grib2_id=(0, 2, 8), units='Pa s**-1', amip='wap', description='The Lagrangian tendency of air pressure, often called "omega", plays the role of the upward component of air velocity when air pressure is being used as the vertical coordinate. If the vertical air velocity is upwards, it is negative when expressed as a tendency of air pressure; downwards is positive. Air pressure is the force per unit area which would be exerted when the moving gas molecules of which the air is composed strike a theoretical surface of any orientation.'), MetVariable(short_name='tau_cirrus', standard_name='tau_cirrus', long_name='Cirrus optical depth', level_type=None, ecmwf_id=None, grib1_id=None, grib2_id=None, units='dimensionless', amip=None, description=None))¶
Minimal set of met variables needed to run the model after pre-processing. The intention here is that
ciwcis unnecessary aftertau_cirrushas already been calculated.
- rad¶
Radiation data formatted as a
MetDatasetat a single pressure level [-1]
- rad_variables = ((MetVariable(short_name='rst', standard_name='toa_net_downward_shortwave_flux', long_name='TOA net downward shortwave flux', level_type='nominalTop', ecmwf_id=None, grib1_id=None, grib2_id=None, units='W m**-2', amip='rst', description='"shortwave" means shortwave radiation. "toa" means top of atmosphere. "Downward" indicates a vector component which is positive when directed downward (negative upward). Net downward radiation is the difference between radiation from above (downwelling) and radiation from below (upwelling). In accordance with common usage in geophysical disciplines, "flux" implies per unit area, called "flux density" in physics.'), MetVariable(short_name='tsr', standard_name='top_net_solar_radiation', long_name='Top of atmosphere net solar (shortwave) radiation', level_type='nominalTop', ecmwf_id=178, grib1_id=None, grib2_id=(0, 4, 1), units='J m**-2', amip=None, description="This parameter is the incoming solar radiation (also known as shortwave radiation) minus the outgoing solar radiation at the top of the atmosphere. It is the amount of radiation passing through a horizontal plane. The incoming solar radiation is the amount received from the Sun. The outgoing solar radiation is the amount reflected and scattered by the Earth's atmosphere and surfaceSee https://www.ecmwf.int/sites/default/files/elibrary/2015/18490-radiation-quantities-ecmwf-model-and-mars.pdf"), MetVariable(short_name='suswrf', standard_name='toa_upward_shortwave_flux', long_name='Top of atmosphere upward shortwave radiation', level_type='nominalTop', ecmwf_id=None, grib1_id=None, grib2_id=(0, 4, 193), units='W m**-2', amip=None, description='This parameter is the outgoing shortwave (solar) radiation at the nominal top of the atmosphere.')), (MetVariable(short_name='rlut', standard_name='toa_outgoing_longwave_flux', long_name='TOA outgoing longwave_flux', level_type='nominalTop', ecmwf_id=None, grib1_id=None, grib2_id=None, units='W m**-2', amip='rlut', description='"longwave" means longwave radiation. "toa" means top of atmosphere. The TOA outgoing longwave flux is the upwelling thermal radiative flux, often called the "outgoing longwave radiation" or "OLR". In accordance with common usage in geophysical disciplines, "flux" implies per unit area, called "flux density" in physics.'), MetVariable(short_name='ttr', standard_name='top_net_thermal_radiation', long_name='Top of atmosphere net thermal (longwave) radiation', level_type='nominalTop', ecmwf_id=179, grib1_id=None, grib2_id=(0, 5, 5), units='J m**-2', amip=None, description='The thermal (also known as terrestrial or longwave) radiation emitted to space at the top of the atmosphere is commonly known as the Outgoing Longwave Radiation (OLR). The top net thermal radiation (this parameter) is equal to the negative of OLR.See https://www.ecmwf.int/sites/default/files/elibrary/2015/18490-radiation-quantities-ecmwf-model-and-mars.pdf'), MetVariable(short_name='sulwrf', standard_name='toa_upward_longwave_flux', long_name='Top of atmosphere upward longwave radiation', level_type='nominalTop', ecmwf_id=None, grib1_id=None, grib2_id=(0, 5, 193), units='W m**-2', amip=None, description='This parameter is the outgoing longwave (thermal) radiation at the nominal top of the atmosphere.')))¶
Required single-level top of atmosphere radiation variables. Variable keys depend on data source (e.g. ECMWF, GFS).
- require_met()¶
Ensure that
metis a MetDataset.- Returns:
MetDataset– Returns reference tomet. This is helpful for type narrowingmetwhen meteorology is required.- Raises:
ValueError – Raises when
metis None.
- require_source_type(type_)¶
Ensure that
sourceistype_.- Returns:
_Source– Returns reference tosource. This is helpful for type narrowingsourceto specific type(s).- Raises:
ValueError – Raises when
sourceis not_type_.
- set_source(source=None)¶
Attach original or copy of input
sourcetosource.- Parameters:
source (
MetDataset | GeoVectorDataset | Flight | Iterable[Flight] | None) – Parametersourcepassed ineval(). If None, an empty MetDataset with coordinates likemetis set tosource.
See also
- set_source_met(optional=False, variable=None)¶
Ensure or interpolate each required
met_variablesonsource.For each variable in
met_variables, checksourcefor data variable with the same name.For
GeoVectorDatasetsources, try to interpolatemetif variable does not exist.For
MetDatasetsources, try to get data frommetif variable does not exist.- Parameters:
optional (
bool, optional) – Includeoptional_met_variablesvariable (
MetVariable | Sequence[MetVariable] | None, optional) – MetVariable to set, frommet_variables. If None, set all variables inmet_variablesandoptional_met_variablesifoptionalis True.
- Raises:
ValueError – Variable does not exist and
sourceis a MetDataset.
- timesteps¶
Array of
numpy.datetime64time steps for contrail evolution