:py:mod:`km3io.tools` ===================== .. py:module:: km3io.tools Module Contents --------------- Classes ~~~~~~~ .. autoapisummary:: km3io.tools.cached_property km3io.tools.TimeConverter Functions ~~~~~~~~~ .. autoapisummary:: km3io.tools.unfold_indices km3io.tools.to_num km3io.tools.unique km3io.tools.uniquecount km3io.tools.get_w2list_param km3io.tools.fitinf km3io.tools.count_nested km3io.tools.get_multiplicity km3io.tools.best_track km3io.tools.mask km3io.tools.has_jmuon km3io.tools.has_jshower km3io.tools.has_aashower km3io.tools.has_dusjshower km3io.tools.best_jmuon km3io.tools.best_jshower km3io.tools.best_aashower km3io.tools.best_dusjshower km3io.tools.is_cc km3io.tools.usr km3io.tools.is_bit_set km3io.tools.is_3dshower km3io.tools.is_mxshower km3io.tools.is_3dmuon km3io.tools.get_w2list_idx km3io.tools.is_nanobeacon km3io.tools.angle km3io.tools.magnitude km3io.tools.theta km3io.tools.theta_separg km3io.tools.phi km3io.tools.phi_separg km3io.tools.zenith km3io.tools.azimuth km3io.tools.angle_between km3io.tools.unit_vector Attributes ~~~~~~~~~~ .. autoapisummary:: km3io.tools.BASKET_CACHE_SIZE km3io.tools.BASKET_CACHE .. py:data:: BASKET_CACHE_SIZE .. py:data:: BASKET_CACHE .. py:class:: cached_property(function) A simple cache decorator for properties. .. !! processed by numpydoc !! .. py:function:: unfold_indices(obj, indices) Unfolds an index chain and returns the corresponding item .. !! processed by numpydoc !! .. py:function:: to_num(value) Convert a value to a numerical one if possible .. !! processed by numpydoc !! .. py:function:: unique(array, dtype=np.int64) Return the unique elements of an array with a given dtype. The performance is better for pre-sorted input arrays. .. !! processed by numpydoc !! .. py:function:: uniquecount(array, dtype=np.int64) Count the number of unique elements in a jagged Awkward1 array. .. !! processed by numpydoc !! .. py:function:: get_w2list_param(events, generator, param) Get all the values of a specific parameter from the w2list in offline neutrino files. :Parameters: **events** : km3io.offline.OfflineBranch events class in offline neutrino files. **generator** : str the name of the software generating neutrinos, it is either 'genhen' or 'gseagen'. **param** : str the name of the parameters found in w2list as defined in the KM3NeT-Dataformat for both genhen and gseagen. :Returns: awkward.Array array of the values of interest. .. !! processed by numpydoc !! .. py:function:: fitinf(fitparam, tracks) Access fit parameters in tracks.fitinf. :Parameters: **fitparam** : int the fit parameter key according to fitparameters defined in KM3NeT-Dataformat (see km3io.definitions.fitparameters). **tracks** : ak.Array or km3io.rootio.Branch reconstructed tracks with .fitinf attribute :Returns: awkward1.Array awkward array of the values of the fit parameter requested. Missing values are set to NaN. .. !! processed by numpydoc !! .. py:function:: count_nested(arr, axis=0) Count elements in a nested awkward Array. :Parameters: **arr** : awkward1.Array Array of data. Example tracks.fitinf or tracks.rec_stages. **axis** : int, optional axis = 0: to count elements in the outmost level of nesting. axis = 1: to count elements in the first level of nesting. axis = 2: to count elements in the second level of nesting. :Returns: awkward1.Array or int counts of elements found in a nested awkward1 Array. .. !! processed by numpydoc !! .. py:function:: get_multiplicity(tracks, rec_stages) Tracks selection based on specific reconstruction stages. Counts how many tracks with the specific reconstructions stages are found per event. :Parameters: **tracks** : km3io.offline.OfflineBranch tracks or a subste of tracks. **rec_stages** : list Reconstruction stages (the ordering is respected) e.g. [1, 2, 3, 4, 5]. :Returns: awkward1.Array tracks multiplicty. .. !! processed by numpydoc !! .. py:function:: best_track(tracks, startend=None, minmax=None, stages=None) Best track selection. :Parameters: **tracks** : awkward.Array A list of tracks or doubly nested tracks, usually from OfflineReader.events.tracks or subarrays of that, containing recunstructed tracks. **startend: tuple(int, int), optional** The required first and last stage in tracks.rec_stages. **minmax: tuple(int, int), optional** The range (minimum and maximum) value of rec_stages to take into account. **stages** : list or set, optional - list: the order of the rec_stages is respected. - set: a subset of required stages; the order is irrelevant. :Returns: awkward.Array or namedtuple Be aware that the dimensions are kept, which means that the final track attributes are nested when multiple events are passed in. If a single event (just a list of tracks) is provided, a named tuple with a single track and flat attributes is created. :Raises: ValueError When invalid inputs are specified. .. !! processed by numpydoc !! .. py:function:: mask(arr, sequence=None, startend=None, minmax=None, atleast=None) Return a boolean mask which mask each nested sub-array for a condition. :Parameters: **arr** : awkward.Array with ndim>=2 The array to mask. **startend: tuple(int, int), optional** True for entries where the first and last element are matching the tuple. **minmax: tuple(int, int), optional** True for entries where each element is within the min-max-range. **sequence** : list(int), optional True for entries which contain the exact same elements (in that specific order) **atleast** : list(int), optional True for entries where at least the provided elements are present. **An extensive discussion about this implementation can be found at:** .. **https://github.com/scikit-hep/awkward-1.0/issues/580** .. **Many thanks for Jim for the fruitful discussion and the final implementation.** .. .. !! processed by numpydoc !! .. py:function:: has_jmuon(tracks) Check if given tracks contain JMUON reconstruction. .. !! processed by numpydoc !! .. py:function:: has_jshower(tracks) Check if given tracks contain JSHOWER reconstruction. .. !! processed by numpydoc !! .. py:function:: has_aashower(tracks) Check if given tracks contain AASHOWER reconstruction. .. !! processed by numpydoc !! .. py:function:: has_dusjshower(tracks) Check if given tracks contain AASHOWER reconstruction. .. !! processed by numpydoc !! .. py:function:: best_jmuon(tracks) Select the best JMUON track. .. !! processed by numpydoc !! .. py:function:: best_jshower(tracks) Select the best JSHOWER track. .. !! processed by numpydoc !! .. py:function:: best_aashower(tracks) Select the best AASHOWER track. .. !! processed by numpydoc !! .. py:function:: best_dusjshower(tracks) Select the best DISJSHOWER track. .. !! processed by numpydoc !! .. py:function:: is_cc(fobj) Determin if events are a result of a Charged Curent interaction (CC) or a Neutral Curent interaction (NC). :Parameters: **fobj** : km3io.offline.OfflineReader offline neutrino file object. :Returns: awkward1.highlevel.Array a mask where True corresponds to an event resulting from a CC interaction. :Raises: ValueError if the simulations program used to generate the neutrino file is neither gseagen nor genhen. .. !! processed by numpydoc !! .. py:function:: usr(objects, field) Return the usr-data for a given field. :Parameters: **objects** : awkward.Array Events, tracks, hits or whatever objects which have usr and usr_names fields (e.g. OfflineReader().events). .. !! processed by numpydoc !! .. py:function:: is_bit_set(value, bit_position) Returns true if a bit at the given position is 1. value: int or array([u]int64) The value to check, can be a single value or an array of values. bit_position: int 0 for the first position, 1 for the second etc. .. !! processed by numpydoc !! .. py:function:: is_3dshower(trigger_mask) Returns True if the trigger mask contains the 3D shower flag. :Parameters: **trigger_mask** : int or array(int) A value or an array of the trigger_mask, either of an event, or a hit. .. !! processed by numpydoc !! .. py:function:: is_mxshower(trigger_mask) Returns True if the trigger mask contains the MX shower flag. :Parameters: **trigger_mask** : int or array(int) A value or an array of the trigger_mask, either of an event, or a hit. .. !! processed by numpydoc !! .. py:function:: is_3dmuon(trigger_mask) Returns True if the trigger mask contains the 3D muon flag. :Parameters: **trigger_mask** : int or array(int) A value or an array of the trigger_mask, either of an event, or a hit. .. !! processed by numpydoc !! .. py:function:: get_w2list_idx(f) Get the correct w2list_idx for the given file, or None if there is none. :Parameters: **f** : km3io.OfflineReader The file. .. !! processed by numpydoc !! .. py:function:: is_nanobeacon(trigger_mask) Returns True if the trigger mask contains the nano-beacon flag. :Parameters: **trigger_mask** : int or array(int) A value or an array of the trigger_mask, either of an event, or a hit. .. !! processed by numpydoc !! .. py:class:: TimeConverter(event) Auxiliary class to convert Monte Carlo hit times to DAQ/triggered hit times. .. !! processed by numpydoc !! .. py:attribute:: FRAME_TIME_NS :value: 100000000.0 .. py:method:: get_time_of_frame(frame_index) Get start time of frame in ns since start of run for a given frame index :Parameters: **frame_index: int** The index of the DAQ frame .. !! processed by numpydoc !! .. py:method:: get_DAQ_time(t0) Get DAQ/triggered hit time :Parameters: **t0: float or array(float)** Simulated time [ns] .. !! processed by numpydoc !! .. py:method:: get_MC_time(t0) Get Monte Carlo hit time :Parameters: **t0: float or array(float)** DAQ/trigger hit time [ns] .. !! processed by numpydoc !! .. py:function:: angle(v1, v2, normalized=False) Compute the unsigned angle between two vectors. For a stacked input, the angle is computed pairwise. Inspired by the "vg" Python package. :Parameters: **v1** : np.array With shape `(3,)` or a `kx3` stack of vectors. **v2** : np.array A vector or stack of vectors with the same shape as `v1`. **normalized** : bool By default, the vectors will be normalised unless `normalized` is `True`. :Returns: A float or a vector of floats with the angle in radians. .. .. !! processed by numpydoc !! .. py:function:: magnitude(v) Calculates the magnitude of a vector or array of vectors. :Parameters: **v** : np.array With shape `(3,)` or `kx3` stack of vectors. :Returns: A float or a vector of floats with the magnitudes. .. .. !! processed by numpydoc !! .. py:function:: theta(v) Neutrino direction in polar coordinates. :Parameters: **v** : array (x, y, z) .. .. rubric:: Notes Downgoing event: theta = 180deg Horizont: 90deg Upgoing: theta = 0 Angles in radians. .. !! processed by numpydoc !! .. py:function:: theta_separg(dir_z) .. py:function:: phi(v) Neutrino direction in polar coordinates. :Parameters: **v** : array (x, y, z) .. .. rubric:: Notes ``phi``, ``theta`` is the opposite of ``zenith``, ``azimuth``. Angles in radians. .. !! processed by numpydoc !! .. py:function:: phi_separg(dir_x, dir_y) .. py:function:: zenith(v) Return the zenith angle in radians. :Parameters: **v** : array (x, y, z) .. .. rubric:: Notes Defined as 'Angle respective to downgoing'. Downgoing event: zenith = 0 Horizont: 90deg Upgoing: zenith = 180deg .. !! processed by numpydoc !! .. py:function:: azimuth(v) Return the azimuth angle in radians. :Parameters: **v** : array (x, y, z) .. .. rubric:: Notes ``phi``, ``theta`` is the opposite of ``zenith``, ``azimuth``. This is the 'normal' azimuth definition -- beware of how you define your coordinates. KM3NeT defines azimuth differently than e.g. SLALIB, astropy, the AAS.org .. !! processed by numpydoc !! .. py:function:: angle_between(v1, v2, axis=0) Returns the angle in radians between vectors 'v1' and 'v2'. If axis=1 it evaluates the angle between arrays of vectors. .. rubric:: Examples >>> angle_between((1, 0, 0), (0, 1, 0)) 1.5707963267948966 >>> angle_between((1, 0, 0), (1, 0, 0)) 0.0 >>> angle_between((1, 0, 0), (-1, 0, 0)) 3.141592653589793 >>> v1 = np.array([[1, 0, 0], [1, 0, 0]]) >>> v2 = np.array([[0, 1, 0], [-1, 0, 0]]) >>> angle_between(v1, v2, axis=1) array([1.57079633, 3.14159265]) .. !! processed by numpydoc !! .. py:function:: unit_vector(vector, **kwargs) Returns the unit vector of the vector. .. !! processed by numpydoc !!