The trajectory_compute module¶
Introduction¶
The objective of this code is to construct sets of back and forward trajectories that pass through specific regions, or objects, at a reference time, and ‘families’ of such sets with different reference times.
The basic class is trajectory_compute.Trajectories
.
However, the class trajectory_compute.Trajectory_Family
is more powerful; it produces a related set, or family, of trajectories with sequential reference times so that objects which interact can be identified. Results can be plotted and animated using The trajectory_plot module.
Defining Objects¶
As of version 0.2, the objects are specified by providing two functions:
- ref_func(dataset, time_index)
Function to return reference trajectory positions and labels. This has two parameters: dataset, a Netcdf file handle containing data at the reference time, and time_index, the index of required time in file. It must return three variables. An array giving positions of trajectory origin points, a corresponding array of point labels (integers) and an integer giving the total number of objects.
- in_obj_func(traj, tr_time, obj_ptrs)
Function to determine which points are inside an object. This is passed a trajectory object and returns a logical array selecting points which are inside objects, and a data array corresponding to the object data criterion (e.g. cloud water). If it is also passed a time index and logical array selecting an object, it these apply to this time and object, otherwise the whole trajectory data array is used.
Both are passed a dictionary of keyword arguments to these functions via kwargs.
An example pair of functions has been provided: trajectory_compute.trajectory_cloud_ref()
and trajectory_compute.in_cloud()
.
New at 0.4
An alternative pair based upon George Efstathiou’s code to study coherent structures in the boundary layer is included. See trajectory_compute.trajectory_cstruct_ref()
and trajectory_compute.in_cstruc()
(which both use trajectory_compute.cstruct_select()
).
An example of their use is:
import glob
import os
import numpy as np
from trajectory_compute import *
dir = 'C:/Users/paclk/OneDrive - University of Reading/traj_data/r16/'
files = glob.glob(dir+"diagnostics_3d_ts_*.nc")
files.sort(key=file_key)
ref_prof_file = glob.glob(dir+'diagnostics_ts_*.nc')[0]
dt = 60
first_ref_time = 50*dt
last_ref_time = 90*dt
tr_back_len = 40*dt
tr_forward_len = 30*dt
ref = 40
kwa={'thresh':1.0E-5}
tfm = trajectory_family(files, ref_prof_file, \
first_ref_time, last_ref_time, \
tr_back_len, tr_forward_len, \
100.0, 100.0, 40.0, \
trajectory_cloud_ref, in_cloud, \
kwargs=kwa )
Using Trajectories objects¶
The most important attributes of the Trajectories class are: trajectory, data, traj_times, labels and nobjects. The trajectory array gives the trajectory positions, and is dimensioned [nt, m, 3] where nt is total number of times, m the number of trajectory points at a given time and last index gives x,y,z (index 0, 1, 2). The data array contains the associated model variables interpolated to the trajectory and is dimensioned [nt, m, n] where n is the number of variables in variable_list.
The variable list can be supplied to the trajectory_family and trajectory classes but default to the following:
variable_list = { \
"u":r"$u$ m s$^{-1}$", \
"v":r"$v$ m s$^{-1}$", \
"w":r"$w$ m s$^{-1}$", \
"th":r"$\theta$ K", \
"p":r"Pa", \
"q_vapour":r"$q_{v}$ kg/kg", \
"q_cloud_liquid_mass":r"$q_{cl}$ kg/kg", \
}
Note this is a dictionary with keys the NetCDF variable names in the MONC output, pointing to format strings for plotting. Unfortunately, the data array is not in the same order as the variable list, due to a quirk in the way python orders dictionaries. The trajectory_compute.Trajectories.var()
method is provided to convert a variable name to an index.
New at 0.4
Some derived variables can be specfied, such as ‘th_L’.
Variables can have have operators (such as spatial derivatives, deviation from mean) applied.
See documentation for trajectory_compute.get_data()
for more detail.
Users are encouraged to expand this functionality, or suggest required expansions.
The traj_times array dimensioned[nt] with times corresponding to trajectory.
The labels array is dimensioned [m] and identifies which object trajectory points and associated data with labels 0 to nobjects-1. Thus, in a trajectory object traj, the points corresponding to object iobj can be found using traj.trajectory[:, traj.labels==iobj, :], and the associated data traj.data[:, traj.labels==iobj, :].
Detailed Module Contents¶
The entire module is documented below.
-
class
trajectory_compute.
Trajectory_Family
(files, ref_prof_file, first_ref_time, last_ref_time, back_len, forward_len, deltax, deltay, deltaz, ref_func, in_obj_func, kwargs={}, variable_list=None)¶ Class defining a family of back trajectories.
This is an ordered list of trajectories with sequential reference times.
- Parameters
files – ordered list of files used to generate trajectories
ref_prof_file – name of file containing reference profile.
first_ref_time – Time of reference.
last_ref_time – Time of reference.
back_len – Time to go back from ref.
forward_len – Time to go forward from ref.
deltax – Model x grid spacing in m.
deltay – Model y grid spacing in m.
deltaz – Model z grid spacing in m.
variable_list=None – List of variable names for data to interpolate to trajectory.
ref_func – function to return reference trajectory positions and labels.
in_obj_func – function to determine which points are inside an object.
kwargs – any additional keyword arguments to ref_func (dict).
-
family
¶ List of trajectory objects with required reference times.
- Type
list
@author: Peter Clark
-
matching_object_list
(master_ref=None, select=None)¶ Method to generate a list of matching objects at all times they match. Matching is done using overlap of in_obj boxes.
- Parameters
master_ref=None – Which trajectory set in family to use to find matching objects. Default is the last set.
select=None – Which objects to find matches for.
- Returns
dictionary pointing to arrays of mean properties and meta data:
Dictionary keys:: "master_ref" : master_ref actually used. "objects" : object numbers selected. "matching_objects" : List described below.
List with one member for each trajectory set with an earlier reference time than master_ref.
By default, this is all the sets apart from the last. Let us say t_off is the backwards offset from the reference time, with t_off=0 for the time immediately before the reference time.
Each member is a list with one member for each back trajectory time in the master_ref trajectories, containing a list with one member for every object in master_ref,comprising an array of integer indices of the matching objects in the earlier set that match the reference object at the given time. This may be empty if no matches have been found. Let us say the it_back measures the time backwards from master_ref, with it_back=0 at the reference time.
Thus, if mol is the matching object list, mol[t_off][it_back][iobj] is an array of object ids belonging to trajectory set master_ref-(t_off+1), at master_ref trajectory time it_back before the master_ref reference time and matching object index iobj a the reference time.
@author: Peter Clark
-
print_matching_object_list
(master_ref=None, select=None)¶ Method to print matching object list.
See method matching_object_list.
@author: Peter Clark
-
matching_object_list_summary
(master_ref=None, select=None, overlap_thresh=0.02)¶ Method to classify matching objects.
- Parameters
master_ref=None – Which trajectory set in family to use to find matching objects. Default is the last set.
Threshold for overlap to be sufficient for inclusion. (overlap_thresh=0.02) –
- Returns
Dictionary pointing to arrays of mean properties and meta data:
Dictionary keys:: "master_ref" : master_ref actually used. "objects" : object numbers selected. "matching_object_summary" : List described below.
List with one member for each trajectory set with an earlier reference than master_ref.
By default, this is all the sets apart from the last.
Each member is a list with one member for each object in master_ref, containing a list of objects in the earlier set that match the object in master_ref AT ANY TIME.
Each of these is classified ‘Linked’ if the object is in the max_at_ref list at the earler reference time, otherwise ‘Same’.
@author: Peter Clark
-
print_matching_object_list_summary
(master_ref=None, select=None, overlap_thresh=0.02)¶ Method to print matching object list summary.
See method matching_object_list_summary.
@author: Peter Clark
-
find_linked_objects
(master_ref=None, select=None, overlap_thresh=0.02)¶ Method to find all objects linked to objects in max_at_ref list in master_ref.
- Parameters
master_ref=None – Which trajectory set in family to use to find matching objects. Default is the last set.
overlap_thresh=0.02 – Threshold for overlap to be sufficient for inclusion.
- Returns
List with one member for each object in max_at_ref list in master_ref. Each member is a an array of triplets containing the time, object id and percentage ovelap with next of objects in the max_at_ref list of the family at time classified as ‘Linked’.
@author: Peter Clark
-
print_linked_objects
(master_ref=None, select=None, overlap_thresh=0.02)¶ Method to print linked object list.
See method find_linked_objects.
@author: Peter Clark
-
find_super_objects
(master_ref=None, overlap_thresh=0.02)¶ Method to find all objects linked contiguously to objects in max_at_ref list in master_ref.
- Parameters
master_ref=None – Which trajectory set in family to use to find matching objects. Default is the last set.
overlap_thresh=0.02 – Threshold for overlap to be sufficient for inclusion.
- Returns
super_objects, len_sup
super_objects is a list with one member for each object in max_at_ref list in master_ref. Each member is an array of triplets containing the time, object id and percentage ovelap with next of objects in the max_at_ref list of the family at time classified as ‘Linked’.
len_sup is an array of the lengths of each super_objects array
@author: Peter Clark
-
refine_object_overlap
(t_off, time, obj, mobj, master_ref=None)¶ - Method to estimate degree of overlap between two trajectory objects.
Reference object is self.family trajectory object obj at time master_ref-time Comparison object is self.family trajectory object mobj at same true time, reference time master_ref-(t_off+1).
- Parameters
t_off (integer) – Reference object is at time index master_ref-time
time (integer) – Comparison object is at time index master_ref-(t_off+1)
obj (integer) – Reference object id
mobj (integer) – Comparison object if.
master_ref=None (integer) – default is last set in family.
- Returns
Fractional overlap
@author: Peter Clark
-
class
trajectory_compute.
Trajectories
(files, ref_prof_file, start_time, ref, end_time, deltax, deltay, deltaz, ref_func, in_obj_func, kwargs={}, variable_list=None)¶ Class defining a set of back trajectories with a given reference time.
This is an ordered list of trajectories with sequential reference times.
- Parameters
files (ordered list of strings.) – Files used to generate trajectories
ref_prof_file (string) – name of file containing reference profile.
start_time (int) – Time for origin of back trajectories.
ref (int) – Reference time of trajectories.
end_time (int) – Time for end of forward trajectories.
deltax (float) – Model x grid spacing in m.
deltay (float) – Model y grid spacing in m.
deltaz (float) – Model z grid spacing in m.
variable_list=None (List of strings.) – Variable names for data to interpolate to trajectory.
ref_func (function) – Return reference trajectory positions and labels.
in_obj_func (function) – Determine which points are inside an object.
kwargs (dict) – Any additional keyword arguments to ref_func.
-
refprof
¶ Reference profile:
'rho' (array) : Reference profile of density. 'p' (array) : Reference profile of pressure. 'th' (array) : Reference profile of potential temperature. 'pi' (array) : Reference profile of Exner pressure.
- Type
Dict
-
data
¶ Data associated with trajectory:
nt is total number of times, m the number of trajectory points at a given time and n is the number of variables in variable_list.
- Type
float array [nt, m, n]
-
trajectory
¶ - Type
Array [nt, m, 3] where the last index gives x,y,z
-
traj_error
¶ - Type
Array [nt, m, 3] with estimated error in trajectory.
-
traj_times
¶ - Type
Array [nt] with times corresponding to trajectory.
-
labels
¶ - Type
Array [m] labelling points with labels 0 to nobjects-1.
-
nobjects
¶ - Type
Number of objects.
-
coords
¶ Model grid info:
'xcoord': grid xcoordinate of model space. 'ycoord': grid ycoordinate of model space. 'zcoord': grid zcoordinate of model space. 'ycoord': grid ycoordinate of model space. 'zcoord': grid zcoordinate of model space. 'deltax': Model x grid spacing in m. 'deltay': Model y grid spacing in m. 'z' : Model z grid m. 'zn' : Model z grid m.
- Type
Dict
-
nx
¶ length of xcoord
- Type
int
-
ny
¶ length of ycoord
- Type
int
-
nz
¶ length of zcoord
- Type
int
-
deltat
¶ time spacing in trajectories.
- Type
float
-
ref_func
¶ - Type
function to return reference trajectory positions and labels.
-
in_obj_func
¶ - Type
function to determine which points are inside an object.
-
ref_func_kwargs
¶ Any additional keyword arguments to ref_func (dict).
- Type
dict
-
files
¶ - Type
Input file list.
-
ref
¶ - Type
Index of reference time in trajectory array.
-
end
¶ - Type
Index of end time in trajectory array. (start is always 0)
-
ntimes
¶ - Type
Number of times in trajectory array.
-
npoints
¶ - Type
Number of points in trajectory array.
-
variable_list
¶ - Type
variable_list corresponding to data.
-
data_mean
¶ - Type
mean of in_obj points data.
-
num_in_obj
¶ - Type
number of in_obj points.
-
centroid
¶ - Type
centroid of in_objy points
-
bounding_box
¶ - Type
box containing all trajectory points.
-
in_obj_box
¶ - Type
box containing all in_obj trajectory points.
-
max_at_ref
¶ - Type
list of objects which reach maximum LWC at reference time.
-
@author
- Type
Peter Clark
-
var
(v)¶ Method to convert variable name to numerical pointer.
- Parameters
v (string) – variable name.
- Returns
Numerical pointer to data in data array.
@author: Peter Clark
-
select_object
(iobj)¶ Method to find trajectory and associated data corresponding to iobj.
- Parameters
iobj (integer) – object id .
- Returns
trajectory_array, associated _data.
@author: Peter Clark
-
trajectory_compute.
dict_to_index
(variable_list, v)¶ Method to convert variable name to numerical pointer.
- Parameters
v (string) – variable name.
- Returns
Numerical pointer to v in variable list.
@author: Peter Clark
-
trajectory_compute.
compute_trajectories
(files, start_time, ref_time, end_time, deltax, deltay, deltaz, variable_list, refprof, ref_func, kwargs={})¶ Function to compute forward and back trajectories plus associated data.
- Parameters
files – Ordered list of netcdf files containing 3D MONC output.
start_time – Time corresponding to end of back trajectory.
ref_time – Time at which reference objects are defined.
end_time – Time corresponding to end of forward trajectory.
variable_list – List of variables to interpolate to trajectory points.
refprof – ref profile.
- Returns
Set of variables defining trajectories:
data_val: Array [nt, m, n] where nt is total number of times, m the number of trajectory points at a given time and n is the number of variables in variable_list. trajectory: Array [nt, m, 3] where the last index gives x,y,z traj_error: Array [nt, m, 3] with estimated error in trajectory. traj_times: Array [nt] with times corresponding to trajectory. labels: Array [m] labelling points with labels 0 to nobjects-1. nobjects: Number of objects. coords : Dict containing model grid info. 'xcoord': grid xcoordinate of model space. 'ycoord': grid ycoordinate of model space. 'zcoord': grid zcoordinate of model space. 'ycoord': grid ycoordinate of model space. 'zcoord': grid zcoordinate of model space. 'deltax': Model x grid spacing in m. 'deltay': Model y grid spacing in m. 'z' : Model z grid m. 'zn' : Model z grid m. deltat: time spacing in trajectories.
@author: Peter Clark
-
trajectory_compute.
extract_pos
(nx, ny, dat)¶ Function to extract 3D position from data array.
- Parameters
nx – Number of points in x direction.
ny – Number of points in y direction.
dat – Array[m,n] where n>=5 if cyclic_xy or 3 if not.
- Returns
Array[m,3] n_pvar : Number of dimensions in input data used for pos.
- Return type
pos
-
trajectory_compute.
trajectory_init
(dataset, time_index, variable_list, deltax, deltay, deltaz, refprof, traj_pos)¶ Function to set up origin of back and forward trajectories.
- Parameters
dataset – Netcdf file handle.
time_index – Index of required time in file.
variable_list – List of variable names.
refprof – Dict with reference theta profile arrays.
traj_pos – array[n,3] of initial 3D positions.
- Returns
Trajectory variables:
trajectory : position of origin point. data_val : associated data so far. traj_error : estimated trajectory errors so far. traj_times : trajectory times so far. coords : Dict containing model grid info. 'xcoord': grid xcoordinate of model space. 'ycoord': grid ycoordinate of model space. 'zcoord': grid zcoordinate of model space. 'deltax': Model x grid spacing in m. 'deltay': Model y grid spacing in m. 'z' : Model z grid m. 'zn' : Model z grid m.
@author: Peter Clark
-
trajectory_compute.
back_trajectory_step
(dataset, time_index, variable_list, refprof, coords, trajectory, data_val, traj_error, traj_times)¶ Function to execute backward timestep of set of trajectories.
- Parameters
dataset – netcdf file handle.
time_index – time index in netcdf file.
variable_list – list of variable names.
refprof – Dict with reference theta profile arrays.
coords – Dict containing 1D arrays giving coordinate spaces of data.
trajectory – trajectories so far. trajectory[0] is position of earliest point.
data_val – associated data so far.
traj_error – estimated trajectory errors to far.
traj_times – trajectory times so far.
- Returns
Inputs updated to new location:
trajectory, data_val, traj_error, traj_times
@author: Peter Clark
-
trajectory_compute.
forward_trajectory_step
(dataset, time_index, variable_list, refprof, coords, trajectory, data_val, traj_error, traj_times, vertical_boundary_option=1)¶ Function to execute forward timestep of set of trajectories.
- Parameters
dataset – netcdf file handle.
time_index – time index in netcdf file.
variable_list – list of variable names.
refprof – Dict with reference theta profile arrays.
coords – Dict containing 1D arrays giving coordinate spaces of data.
trajectory – trajectories so far. trajectory[-1] is position of latest point.
data_val – associated data so far.
traj_error – estimated trajectory errors to far.
traj_times – trajectory times so far.
- Returns
Inputs updated to new location:
trajectory, data_val, traj_error, traj_times
@author: Peter Clark
-
trajectory_compute.
whichbox
(xvec, x)¶ Find ix such that xvec[ix]<x<=xvec[ix+1]. Uses numpy.searchsorted.
- Parameters
xvec – Ordered array.
x – Values to locate in xvec.
- Returns
Array of indices.
@author : Peter Clark
-
trajectory_compute.
tri_lin_interp
(data, varp_list, pos, xcoord, ycoord, zcoord, maxindex=None)¶ Tri-linear interpolation with cyclic wrapround in x and y.
- Parameters
data – list of Input data arrays on 3D grid.
varp_list – list of grid info lists.
pos (Array[n,3]) – Positions to interpolate to (in grid units) .
xcoord – 1D coordinate vector.
ycoord – 1D coordinate vector.
zcoord – 1D coordinate vector.
- Returns
list of 1D arrays of data interpolated to pos.
@author : Peter Clark
-
trajectory_compute.
data_to_pos
(data, varp_list, pos, xcoord, ycoord, zcoord, maxindex=None)¶ Function to interpolate data to pos.
- Parameters
data – list of data array.
varp_list – list of grid info lists.
pos – array[n,3] of n 3D positions.
xcoord,ycoord,zcoord – 1D arrays giving coordinate spaces of data.
- Returns
list of arrays containing interpolated data.
@author: Peter Clark
-
trajectory_compute.
load_traj_pos_data
(dataset, it)¶ Function to read trajectory position variables from file. :param dataset: netcdf file handle. :param it: time index in netcdf file.
- Returns
List of arrays containing interpolated data.
@author: Peter Clark
-
trajectory_compute.
d_by_dx_field
(field, dx, varp, xaxis=0)¶ Numerically differentiate field in x direction assuming cyclic, mofiying grid descriptor.
- Parameters
field (array) – Data array, usualy 3D (x, y, z) but should be general.
dx (float) – Grid spacing
varp (List of logicals) – Describes grid offset. True means offset from p on grid in that axis.
xaxis (int) – Pointer to x dimension - usually 0 but can be set to, e.g., 1 if time- dimension included in data.
- Returns
- Return type
None.
-
trajectory_compute.
d_by_dy_field
(field, dy, varp, yaxis=1)¶ Numerically differentiate field in y direction assuming cyclic, mofiying grid descriptor.
- Parameters
field (array) – Data array, usualy 3D (x, y, z) but should be general.
dx (float) – Grid spacing
varp (List of logicals) – Describes grid offset. True means offset from p on grid in that axis.
xaxis (int) – Pointer to y dimension - usually 1 but can be set to, e.g., 2 if time- dimension included in data.
- Returns
- Return type
None.
-
trajectory_compute.
d_by_dz_field
(field, z, zn, varp)¶ Numerically differentiate field in z direction, mofiying grid descriptor.
- Parameters
field (array) – Data array, usualy 3D (x, y, z) but should be general. Assumes z axis is last dimension.
dx (float) – Grid spacing
varp (List of logicals) – Describes grid offset. True means offset from p on grid in that axis.
xaxis (int) – Pointer to y dimension - usually 1 but can be set to, e.g., 2 if time- dimension included in data.
- Returns
- Return type
None.
-
trajectory_compute.
padleft
(f, zt, axis=0)¶ Add dummy field at bottom of nD array
- Parameters
f – nD field
zt – 1D zcoordinates
axis=0 – Specify axis to extend
- Returns
extended field, extended coord
@author: Peter Clark
-
trajectory_compute.
padright
(f, zt, axis=0)¶ Add dummy field at top of nD array
- Parameters
f – nD field
zt – 1D zcoordinates
axis=0 – Specify axis to extend
- Returns
extended field, extended coord
@author: Peter Clark
-
trajectory_compute.
get_data
(source_dataset, var_name, it, refprof, coords)¶ Extract data from source NetCDF dataset, derived and/or processed data.
Currently supported derived data are:
'th_L' : Liquid water potential temperature. 'th_v' : Virtual potential temperature. 'q_total' : Total water 'buoyancy': Bouyancy based on layer-mean theta_v.
Currently supported operators are:
'v_prime' : Deviation from level mean for any variable v. 'v_crit' : Level critical value for any variable v. 'dv_dx' : Derivative of v in x direction. 'dv_dy' : Derivative of v in y direction. 'dv_dz' : Derivative of v in z direction.
Operators must be combined using parentheseses, e.g. ‘d(d(th_prime)_dx)_dy’.
- Parameters
source_dataset – handle of NetCDF dataset.
var_name – String describing data required.
it – Time index required.
refprof – Dict containg reference profile.
coords – Dict containing model coords.
- Returns
variable, variable_grid_properties
@author: Peter Clark and George Efstathiou
-
trajectory_compute.
load_traj_step_data
(dataset, it, variable_list, refprof, coords)¶ Function to read trajectory variables and additional data from file for interpolation to trajectory.
- Parameters
dataset – netcdf file handle.
it – time index in netcdf file.
variable_list – List of variable names.
refprof – Dict with reference theta profile arrays.
- Returns
List of arrays containing interpolated data.
@author: Peter Clark
-
trajectory_compute.
phase
(vr, vi, n)¶ Function to convert real and imaginary points to location on grid size n.
- Parameters
vr,vi – real and imaginary parts of complex location.
n – grid size
- Returns
Real position in [0,n)
@author: Peter Clark
-
trajectory_compute.
label_3D_cyclic
(mask)¶ Function to label 3D objects taking account of cyclic boundary in x and y. Uses ndimage(label) as primary engine.
- Parameters
mask – 3D logical array with object mask (i.e. objects are contiguous True).
- Returns
Object identifiers:
labs : Integer array[nx,ny,nz] of labels. -1 denotes unlabelled. nobjs : number of distinct objects. Labels range from 0 to nobjs-1.
@author: Peter Clark
-
trajectory_compute.
unsplit_object
(pos, nx, ny)¶ - Function to gather together points in object separated by cyclic boundaries.
For example, if an object spans the 0/nx boundary, so some points are close to zero, some close to nx, they will be adjusted to either go from negative to positive, close to 0, or less than nx to greater than. The algorithm tries to group on the larges initial set. Uses sklearn.cluster.KMeans as main engine.
- Parameters
pos – grid positions of points in object.
nx,ny – number of grid points in x and y directions.
- Returns
Adjusted grid positions of points in object.
@author: Peter Clark
-
trajectory_compute.
unsplit_objects
(trajectory, labels, nobjects, nx, ny)¶ Function to unsplit a set of objects at a set of times using unsplit_object on each.
- Parameters
trajectory – Array[nt, np, 3] of trajectory points, with nt times and np points.
labels – labels of trajectory points.
nx,ny – number of grid points in x and y directions.
- Returns
Trajectory array with modified positions.
@author: Peter Clark
-
trajectory_compute.
compute_traj_boxes
(traj, in_obj_func, kwargs={})¶ Function to compute two rectangular boxes containing all and in_obj points for each trajectory object and time, plus some associated data.
- Parameters
traj – trajectory object.
in_obj_func – function to determine which points are inside an object.
kwargs – any additional keyword arguments to ref_func (dict).
- Returns
Properties of in_obj boxes:
data_mean : mean of points (for each data variable) in in_obj. in_obj_data_mean : mean of points (for each data variable) in in_obj. num_in_obj : number of in_obj points. traj_centroid : centroid of in_obj box. in_obj__centroid : centroid of in_obj box. traj_box : box coordinates of each trajectory set. in_obj_box : box coordinates for in_obj points in each trajectory set.
@author: Peter Clark
-
trajectory_compute.
print_boxes
(traj)¶ Function to print information provided by compute_traj_boxes stored as attributes of traj.
- Parameters
traj – trajectory object.
@author: Peter Clark
-
trajectory_compute.
box_overlap_with_wrap
(b_test, b_set, nx, ny)¶ Function to compute whether rectangular boxes intersect.
- Parameters
b_test – box for testing array[8,3]
b_set – set of boxes array[n,8,3]
nx – number of points in x grid.
ny – number of points in y grid.
- Returns
indices of overlapping boxes
@author: Peter Clark
-
trajectory_compute.
find_time_in_files
(files, ref_time, nodt=False)¶ - Function to find file containing data at required time.
Assumes file names are of form *_tt.0* where tt is model output time.
- Parameters
files – ordered list of files
ref_time – required time.
nodt – if True do not look for next time to get delta_t
- Returns
Variables defining location of data in file list:
ref_file: Index of file containing required time in files. it: Index of time in dataset. delta_t: Interval between data.
@author: Peter Clark
-
trajectory_compute.
trajectory_cloud_ref
(dataset, time_index, thresh=1e-05)¶ Function to set up origin of back and forward trajectories.
- Parameters
dataset – Netcdf file handle.
time_index – Index of required time in file.
thresh=0.00001 – Cloud liquid water threshold for clouds.
- Returns
Trajectory variables:
traj_pos : position of origin point. labels : array of point labels. nobjects : number of objects.
@author: Peter Clark
-
trajectory_compute.
in_cloud
(traj, *argv, thresh=1e-05)¶ Function to select trajectory points inside cloud. In general, ‘in-object’ identifier should return a logical mask and a quantitative array indicating an ‘in-object’ scalar for use in deciding reference times for the object.
- Parameters
traj – Trajectory object.
*argv –
**kwargs –
- Returns
Variables defining those points in cloud:
mask : Logical array like trajectory array. qcl : Cloud liquid water content array.
@author: Peter Clark
-
trajectory_compute.
cloud_select
(qcl, thresh=1e-05)¶ Simple function to select in-cloud data; used in in_cloud and trajectory_cloud_ref for consistency. Written as a function to set structure for more complex object identifiers.
- Parameters
qcl – Cloud liquid water content array.
thresh=0.00001 – Cloud liquid water threshold for clouds.
- Returns
Logical array like trajectory array selecting those points inside cloud.
@author: Peter Clark
-
trajectory_compute.
trajectory_cstruct_ref
(dataset, time_index, thresh=1e-05, find_objects=True)¶ Function to set up origin of back and forward trajectories.
- Parameters
dataset – Netcdf file handle.
time_index – Index of required time in file.
thresh=0.00001 – Cloud liquid water threshold for clouds.
find_objects=True – Enable idenfificatoin of distinct objects.
- Returns
Trajectory variables:
traj_pos : position of origin point. labels : array of point labels. nobjects : number of objects.
@author: George Efstathiou and Peter Clark
-
trajectory_compute.
in_cstruc
(traj, *argv, **kwargs)¶ Function to select trajectory points inside coherent structure. In general, ‘in-object’ identifier should return a logical mask and a quantitative array indicating an ‘in-object’ scalar for use in deciding reference times for the object.
- Parameters
traj – Trajectory object.
*argv –
**kwargs –
- Returns
Variables defining those points in cloud:
mask : Logical array like trajectory array. qcl : Cloud liquid water content array.
@author: George Efstathiou and Peter Clark
-
trajectory_compute.
cstruct_select
(qcl, w, zpos, tr1, tr2, thresh=1e-05, ver=4)¶ Simple function to select in-coherent structure data; used in in_cstruct and trajectory_cstruct_ref for consistency. Written as a function to set structure for more complex object identifiers.
- Parameters
qcl – Cloud liquid water content array.
w – Vertical velocity array.
zpos – Height array.
tr1 – tracer 1 array.
tr2 – tracer 2 array.
thresh=0.00001 – Cloud liquid water threshold for clouds.
- Returns
Logical array like trajectory array selecting those points inside cloud.
@author: George Efstathiou and Peter Clark