waveformtools.waveformtools
This module contains a set of tools to handle data from numerical relativity simulations and gravitational waves.
Functions
|
Append zeros to an array without tapering. |
|
Identify the Approximate start and endpoints of the data. |
|
Function to bin the data and interpolate it at specified width and order. |
|
Center a waveform (wvp, wvc) at the peak of the complex modulous sqrt(wvp**2 + wvc**2). |
|
Check the data (time,datar,datai) for repetetive rows and remove them. |
|
Old version. |
|
Coalign two timeseries. |
|
Coalign two waveforms function 2. |
|
Compute the chirpmass from a2, the coefficient of time of the Finn-Chernoff waform model. |
|
Compute the Newtonian instantaneous frequency of strain waveform from coalescence time and chirp mass, from the Finn-Chernoff model. |
|
Differentiate a timeseries in time domain using the Simple Euler method |
|
Flatten a (3rd order) list of list of lists. |
|
Get the time axis of the waveform centered |
|
Get the approximate starting frequency of the input data by averaging over the first npoints number of points. |
|
Get the angular frequency of the waveform given the complex waveform time step. |
|
Integrate a timeseries using first order method. |
|
Wrapper function for interpolation and resampling. |
|
Function to interpolate a list of timeseries data using the user specified interp_func function and the keyword arguments. |
|
Check if the data has discontinuities. |
|
Check if the data has discontinuities. |
|
Equalize the length of two timeseries/array by appending zeros at the end of the array. |
|
A function to load python objects from the disk using pickle. |
|
Apply low frequency cut filter using a butterworth filter. |
|
Compute the mass ratio from chirpmass. |
|
Match two waveforms and return the time shift, phase shift, normalized waveforms and match coefficient. |
|
Match two waveforms using pycbc subroutines and return the time shift, phase shift, normalized waveforms and match coefficient. |
|
Function to smoothen data. |
|
The print function with verbosity levels and logging facility. |
|
Find the mode of a list |
|
Calculate the norm of a vector. |
|
Calcuate the overlap between two data vectors weighted by the given psd. |
|
A Basic plotting function. |
|
Match function for post merger waveforms. |
|
Display the progress bar to std out from present_count and total_count. |
|
Remove Nans from (xdata,ydata) data pair. |
|
Remove zeros from the input waveform from either sides. |
|
Function to generate timeseries out of the given interpolated data function, epoch,sampling frequency, length(duration). |
|
Resample the waveform pairs. |
|
Roll the data circularly. |
|
A function to save python objects to disk using pickle. |
|
Timeshift an array. |
|
Shorten an array given the start and end points. |
|
Simple match the given waveforms. |
|
Use the Savitzky-Golay Filter to smoothen the data. |
|
Identify the start and endpoints of the data. |
|
A method to taper and append additional zeros at either ends, using the taper function of the pycbc TimeSeries object. |
|
Taper a waveform with a \(tanh\) function at either ends |
|
Taper and equalize the lengths of two arrays. |
|
Find total mass from mass ratio and chirpmass. |
|
Unwrap the phase by finding turning points in phi0. |
|
Given real and imaginary parts of a complex timeseries, extract the amplitude of the complex data vector : (tsdata_p + i * tsdata_x) |
|
Wrapper for extracting the amplitude and the phase of the complex vector. |
|
Given real and imaginary parts of a complex timeseries, extract the phase of the waveform :arctan_(Img(data)/Re(data)) |
- waveformtools.waveformtools.addzeros(data, zeros)[source]
Append zeros to an array without tapering.
- Parameters:
- data1d array or a pycbc TimeSeries
The waveform data as list or numpy array or pycbc timeseries.
- zerosint
The number of zeros to be added.
- Returns:
- data with `zeros` number of zeros concatenated at the end as numpy 1d arraytrial
- waveformtools.waveformtools.apxstartend(data, tol=1e-05)[source]
Identify the Approximate start and endpoints of the data.
- Parameters:
- data: 1d array or a pycbc TimeSeries object
The input waveform whose start and end points need to be identified based on a given tolerance value.
- tol: float
The tolerance value below which the signal is assumed to be absent.
- Returns:
- loc_pair: int (2)
The pair of indices denoting the start and end points of an array
Notes
The starting and ending index of the peak tol (default 1e-5) part of the data is the identification criterion. Requires the data to fall off to tol*peak absolute value outside a certain range.
- waveformtools.waveformtools.bintp(xdata, func_x, width, order, to_plot=True)[source]
Function to bin the data and interpolate it at specified width and order.
- Parameters:
- xdata: 1d array
1D list or numpy array.
- func_x: 1d array
The y axis.
- width: int
Window size for smoothening,
- order: int
Order of the polynomial used for interpolation.
- to_plot: bool
True or False. To plot or not plot the results.
- Returns:
- hist: a list
[binloc, yvals], The location of the bins and the y values associated with the bins.
- waveformtools.waveformtools.center(wvp, wvc=None, delta_t=None)[source]
Center a waveform (wvp, wvc) at the peak of the complex modulous sqrt(wvp**2 + wvc**2).
- Parameters:
- wvp: 1d array or a pycbc TimeSeries object
- wvc: 1d array or a pycbc TimeSeries object
The one/two components of the waveforms as 1d list or numpy arrays or pycbc timeseries.
- delta_t: float
The timestepping delta_t.
- Returns:
- centered_wf: a pycbc TimeSeries objet
The two 1d centered waveform(s) as individual pycbc timeseries.
- waveformtools.waveformtools.cleandata(data, toldt=0.001, bridge='no')[source]
Check the data (time,datar,datai) for repetetive rows and remove them.
- Parameters:
- data: list
Input as a list of 1d arrays [time, data1, data2, …]. All the data share the common time axis time
- toldtfloat
The tolerance for error in checking. defaluts to toldt=1e-3.
- bridgebool
A bridge flag to decide whether or not to interpolate and resample to fill in jump discontinuities.
- Returns:
- cleaned_data: list
The cleaned data array with repetitive rows and gaps (if bridge=True) removed.
- waveformtools.waveformtools.cleandata_old(data, verbose=False)[source]
Old version. Check the data (time,datar,datai) for repetetive rows and remove them.
- Parameters:
- data: list
Input as a list of 1d arrays [time, data1, data2, …]. All the data share the common time axis time
- verbosebool
A verbosity flag.
- Returns:
- cleaned_data: nd array
The cleaned data array with repetitive rows and gaps (if bridge=True) removed
- waveformtools.waveformtools.coalignwfs(tsdata1, tsdata2, delta_t=None)[source]
- Coalign two timeseries.
Wrapper and modification around pycbc functions with some additional functionalities.
- Parameters:
- tsdata1: list, 1d array or a pycbc TimeSeries
- tsdata2: list, 1d array or a pycbc TimeSeries
The two data vectors as 1d lists or numpy arrays or pycbc TimeSeries
- delta_t: float
The time stepping.
- Returns:
- ctsdata1: a pycbc TimeSeries
- tsdata2: a pycbc TimeSeries
A pair of pycbc TimeSeries objects; the aligned first waveform and the second.
- waveformtools.waveformtools.coalignwfs2(tsdata1, tsdata2, delta_t=None)[source]
Coalign two waveforms function 2.
- Parameters:
- tsdata1: 1d array or a pycbc TimeSeries
- tsdata2: 1d array or a pycbc TimeSeries
two data vectors as 1d lists or numpy arrays or pycbc timeseries.
- Returns:
- aligned_waveforms: list
The aligned waveforms in the format [aligned_wf1, aligned_wf2, [norm1, norm2, location of maximum]].
Notes
Procedure:
Adjust length of either waveforms if needed
Normalize the two timeseries
Compute The complex SNR
Shift and roll the first
Returns normalized waveforms
See coalignwfs for description. This algorithm does not use the builtin coalign function from pycbc. This involves normalization of the data vectors explicitly and identifiies the timeshift by computing the complex SNR and finding the maximum element.
- waveformtools.waveformtools.compute_chirp_mass(a2_param)[source]
Compute the chirpmass from a2, the coefficient of time of the Finn-Chernoff waform model.
- Parameters:
- a2_param: 1d array float
The a2 parameter in the Finn Chernoff model.
- Returns:
- chirp_mass: float
The chirp mass
- waveformtools.waveformtools.compute_frequencies(t_coal, t_val, chirp_mass)[source]
Compute the Newtonian instantaneous frequency of strain waveform from coalescence time and chirp mass, from the Finn-Chernoff model.
- Parameters:
- t_coal: float
The coalescence time,
- t_val: 1d array
The time,
- chirp_mass: float
The chirpmass.
- Returns:
- freqs: float
The instantaneous frequency of the strain waveform.
- waveformtools.waveformtools.differentiate(data, delta_t=None, TS=False)[source]
Differentiate a timeseries in time domain using the Simple Euler method
- Parameters:
- data: 1d array a pycbc TimeSeries object
The time series to be differentiated.
- delta_t: float
The grid spacing delta_t. Supplying delta_t overrides the delta_t attribute of TimeSeries.
- Returns:
- ddt_data: pycbc TimeSeries object
The differentiated 1d data as pycbc TimeSeries
- waveformtools.waveformtools.flatten(nflist)[source]
Flatten a (3rd order) list of list of lists. i.e. a three tier list [[[],[]], [[],[]] —> []. This is useful e.g. when combining the data from the list output of multiple MPI ranks.
- Parameters:
- nflist: a third tier list
A list of list of lists (a list of depth three).
- Returns:
- flattened_list: a list
The flattened list i.e. a tier one list.
- waveformtools.waveformtools.get_centered_taxis(time_ax, amps)[source]
- Get the time axis of the waveform centered
at its maximum absolute amplitude.
- Parameters:
- time_ax: 1d array
The 1d array containg the original (uncentered)time axis of the wveform.
- amps1d array
The amplitude of the waveform.
- Returns:
- time_centered: 1d array
The centered time axis of the waveform. The maximum amplitude timestap will be at \(t=0\).
- waveformtools.waveformtools.get_starting_angular_frequency(waveform, delta_t, npoints=400)[source]
Get the approximate starting frequency of the input data by averaging over the first npoints number of points.
- Parameters:
- waveform: 1d array
The 1d complex array of the input waveform.
- delta_t: float
The time step.
- npoints: int
The number of points to average over.
- Returns:
- omega0: float
The approximate starting angular frequency.
Notes
Please suppy a conditioned input waveform that is neatly clipped, and not tapered.
- waveformtools.waveformtools.get_waveform_angular_frequency(waveform, delta_t, timeaxis=None, method='FD')[source]
Get the angular frequency of the waveform given the complex waveform time step. The phase is extracted and is differentiated using one of the available methods to compute the angular frequency.
- Parameters:
- waveform: 1d array
The 1d complex array of the input waveform.
- delta_t: float
The time step.
- timeaxis: 1d array, optional
The time axis of the waveform. Recommended especially when the sampling is non-uniform and Chebyshev method is chosen.
- method: str
The method for computing the derivative. Can be FD or CS. See below for more information.
- Returns:
- omega: 1d array
The real instantaneous frequency of the waveform.
Notes
Available methods
- Chebyshev series (CS): The phase is expanded in
a Chebyshev series and then the differentiated.
- Finite Differencing (FD)A 11 point finite difference
scheme is used to differentiate the phase, and is then smoothened using the Savgol filter.
- waveformtools.waveformtools.integrate_first_order(data, t_start=None, t_end=None, delta_t=None, to_taper=False, TS=False)[source]
Integrate a timeseries using first order method.
- Parameters:
- data: 1d array or a pycbc TimeSeries object
The input function to be integrated.
- delta_t: float
The grid spacing delta_t,
- to_taper: bool
True or False. Whether or not to taper the input series.
- Returns:
- int_data: a pycbc TimeSeries object
TimeSeries of the time integrated data
Notes
Capabilities
Simple Euler integrator, option to taper the ends.
- waveformtools.waveformtools.interp_resam_wfs(wavf_data, old_taxis, new_taxis, resam_kind='cubic')[source]
Wrapper function for interpolation and resampling.
- Parameters:
- wavf_data: 1d array
The yaxis to be interpolated,
- old_taxis, new_taxis: 1darray
Old and New time axis.
- Returns:
- resam_wf_data: 1d array
Interpolated and resampled data.
- waveformtools.waveformtools.interpolate_wfs(ts_data, interp_func, delta_t=None, **kwargs)[source]
Function to interpolate a list of timeseries data using the user specified interp_func function and the keyword arguments.
- Parameters:
- ts_data: list
The 1d data. A list of waveforms as a list or numpy array or pycbc TimeSeries.
- interp_fun: function
An interpolating function.
- delta_t: float
Timestep.
- ``**kwargs``: keyword arguments
additional arguments to the user specified interp_func.
- Returns:
- interp_data: list
A list containing interpolated data.
- waveformtools.waveformtools.iscontinuous(data, delta_t=0, toldt=0.001)[source]
Check if the data has discontinuities. This checks for repetitive time rows and jumps.
- Parameters:
- data: list
Input as a list of 1d arrays [time, data1, data2, …]. All the data share the common time axis time
- delta_t: float
The time stepping.
- toldtfloat
The tolerance for error in checking. defaluts to toldt=1e-3.
- Returns:
- discontinuity_details: a list.
It contains: 1. A list. details of discontinuity: index location of original array,
the corresponding discinbtinuity type.
A float. the global discontinuity type.
Notes
Types of discontunuities
0: Continuous. 1: Repetitive rows. 2: Jumps in timeaxis.
- waveformtools.waveformtools.iscontinuous_old(timeaxis, delta_t=0)[source]
Check if the data has discontinuities. This checks for repetitive time rows and jumps.
- Parameters:
- timeaxis: list
Input as a single 1d time axis or a list of 1d arrays [time, data1, data2, …]. All the data share the common time axis time
- delta_t: float
The time stepping.
- toldtfloat
The tolerance for error in checking. Defaluts to toldt=1e-3.
- Returns:
- discontinuity_detailsa list.
It contains: [ the actual location of discontinuity along the time axis, value of time location of original array, the type of discontinuity].
Notes
Types of discontunuities
0: Continuous. 1: Repetitive rows. 2: Jumps in timeaxis.
- waveformtools.waveformtools.lengtheq(data_a, data_b, delta_t=None, is_ts=False)[source]
Equalize the length of two timeseries/array by appending zeros at the end of the array. No tapering.
- Parameters:
- data_a: list
The input waveform A
- data_b: list
The input waveform B.
- delta_t: float
The time steping. Defaults to None.
- is_ts: bool
To determine whether the given data is a pycbc TimeSeries.
- Returns
- ——-
- equalized_signalslist
The Tapered, length equalized waveforms data_a and data_b, and a flag denoting which waveform was changed, a or b.
Notes
Recommended usage:
Change length of waveform a to match with that of waveform b.
- waveformtools.waveformtools.load_obj(name, obj_dir='./')[source]
A function to load python objects from the disk using pickle.
- Parameters:
- name: string
The filename.
- obj_dirstring
The path to directory in which file exists. Defaults to PWD.
- Returns:
- obj: object
A python object with the contents of the file.
- waveformtools.waveformtools.low_cut_filter(utilde, freqs, order=2, omega0=0.03)[source]
Apply low frequency cut filter using a butterworth filter.
- Parameters:
- utilde: 1d array
The frequency domain data.
- freqs: 1d array
The frequencies.
- order: int
The order of the butterworth filter.
- omega0: float
The cutoff frequency of the butterworth filter.
- Returns:
- utilde_lc: 1d array
The filtered data.
- waveformtools.waveformtools.massratio(chirp_mass)[source]
Compute the mass ratio from chirpmass. Assumes total mass to be 1.
- Parameters:
- mchirp: float
The chirp mass of the system.
- Returns:
- mass_ratio: float
The Mass ratio of the system
- waveformtools.waveformtools.match_wfs(all_time_axes, all_waveforms, delta_t='auto')[source]
Match two waveforms and return the time shift, phase shift, normalized waveforms and match coefficient.
- Parameters:
- time_axes: list
A list containing the time axes of the two waveforms
- waveforms: list
A list of two waveforms. Each is a 1d array.
- delta_t: float, string, optional
The time step of the resampled arrays. Can be A, B, auto or any float value.
- Returns:
- match_details: dict
A dictionary containing the i). match coeffient ii). time_shift iii). phase shift in radians iv). normalized, resampled, waveforms and their
time-axes.
- waveformtools.waveformtools.match_wfs_pycbc(all_time_axes, all_waveforms)[source]
Match two waveforms using pycbc subroutines and return the time shift, phase shift, normalized waveforms and match coefficient.
- Parameters:
- time_axes: list
A list containing the time axes of the two waveforms
- waveforms: list
A list of two waveforms. Each is a 1d array.
- change: int
Which waveform to change, 1 or 2.
- Returns:
- match_details: dict
A dictionary containing the i). match coeffient ii). time_shift iii). phase shift iv). normalized waveforms and their
time-axes.
- waveformtools.waveformtools.mavg(func_x, width)[source]
Function to smoothen data. Moving average over the window width.
- Parameters:
- func_x: 1d array
A list or numpy array of y axis.
- Width: int
The width of the moving average window.
- Returns:
- func_x_avgd: 1d array
1D array of moving averaged y axis.
- waveformtools.waveformtools.message(*args, message_verbosity=2, print_verbosity=2, log_verbosity=2, **kwargs)[source]
The print function with verbosity levels and logging facility.
- Parameters:
- ``*args`non-keyword arguments
same arguments as to that of the print functions,
- message_verbosityint
- print_verbosityint
- log_verbosityint
- ``**kwargs``keyword arguments
Same as that of the print function.
- Returns:
- 1int
messages to stdout and logging of messages, while the function returns 1.
Notes
Verbosity choices:
- message_verbosityEach message carries with it a verbosity level.
More the verbosity less the priority. Default value is 2 i.e. informational.
print_verbosity : prints all messages above this level of verbosity.
log_verbosity : logs all messages above this level of verbosity.
Verbosity levels:
0: Errors
1: Warnings
2: Information
- waveformtools.waveformtools.mode(a_list)[source]
Find the mode of a list
- Parameters:
- a_list: list
The list whose mode is to be found
- Returns:
- list_mode: float
The mode of the list.
- waveformtools.waveformtools.norm(hdat, psd=1.0)[source]
Calculate the norm of a vector.
- Parameters:
- hdat: 1d array or a pycbc TimeSeries object.
The input waveform.
- psa: 1d array
The noise power spectral density of the inner product.
- Returns:
- norm_f: float
The norm with weighting by the psd.
- waveformtools.waveformtools.olap(data1, data2, psd=1)[source]
Calcuate the overlap between two data vectors weighted by the given psd.
- Parameters:
- data1: 1d array or a pycbc TimeSeries object
- data21d array or a pycbc RimeSeries object
The input waveforms
- psd: 1d array
The power spectral density to weight.
- Returns:
- overlap: float
The overlap divided by the psd.
- waveformtools.waveformtools.plot(xdata, func_x, save='no')[source]
A Basic plotting function.
- Parameters:
- xdata: 1d array
The x axis of the function,
- func_x: 1d array
The y axis of the function f(x),
- save: bool.
True or False. Whether the plot should be saved or not.
- Returns:
- 1: (int)
- plots: figures to stdout and disk
Displays the plot, and Saves with the filename provided.
- waveformtools.waveformtools.pmmatch_wfs(waveforms, offset=25, crop=None)[source]
Match function for post merger waveforms.
- Parameters:
- waveforms: a list of pairs
The pairs of waveforms to match in the format [[wf1_pair1, wf2_pair2], [wf1_pair_2, wf2_pair2], …].
- offset: int
Number of indices to shift the data.
- crop: string
A string to decide how to crop the waveforms. The available Options are 1. signal 2. template 3. both.
- Returns:
- matchdet: a list of dicts
A list of dictionaries. Each contains 1. the waveform pair, 2. the match score, 3. the shift index. to maximize the match.
- waveformtools.waveformtools.progressbar(present_count, total_counts, normalize='yes')[source]
Display the progress bar to std out from present_count and total_count.
- Parameters:
- present_count: int
The present count state.
- total_counts: int
The final state.
- Returns:
- 1int
The progress bar is messageed to stdout.
- waveformtools.waveformtools.removeNans(xdata, ydata)[source]
Remove Nans from (xdata,ydata) data pair. Removes Nans in xdata and ydata and the corresponding y and x entries.
- Parameters:
- xdata: 1d array
The x axis of the data.
- ydata1d array
The y axis of the data.
- Returns:
- x_no_nan: 1d array
- y_no_nan: 1d,array
The data pair x,y with Nans removed.
- waveformtools.waveformtools.removezeros(data, delta_t)[source]
Remove zeros from the input waveform from either sides. Similar to startend but return the truncated array.
- Parameters:
- data: 1d array or a pycbc TimeSeries
The input waveform.
- delta_t: float, optional
The time stepping.
- Returns:
- short_ts: a list
A list containing waveforms with zeros removed on either sides, the start and end indices in the format [short_ts, [start_index, end_index]]
- waveformtools.waveformtools.resample(interp_data, new_delta_t, epoch, length, old_delta_t=None)[source]
Function to generate timeseries out of the given interpolated data function, epoch,sampling frequency, length(duration).
- Parameters:
- interp_data: 1d array
The yaxis to be interpolated.
- epoch: float
The starting point in time.
- delta_t: float
New grid spacing to be sampled at.
- length: int
The duration of x axis.
- Returns:
- data: list
A list containing resampled data as pycbc TimeSeries.
- waveformtools.waveformtools.resample_wfs(both_time_axes, both_waveforms, delta_t='auto', Plot=False)[source]
Resample the waveform pairs.
- Parameters:
- both_time_axeslist
A list containing two 1d arrays representing the time axes.
- both_waveformslist
A list containing two 1d arrays respresenting the waveforms.
- delta_tstring, float
The time step to resample at. Auto uses the finest of the two available. A float value can be provided by the user as well.
- Returns:
- both_time_axes_resam1d array
A 1d array representing the resampled time axes.
- both_waveforms_resamlist
A list containing two 1d arrays respresenting the resampled waveforms.
- waveformtools.waveformtools.roll(tsdata, i_roll, is_ts=False)[source]
Roll the data circularly. Circular counterpart of shiftmatched function.
- Parameters:
- tsdata: 1d array or pycbc TimeSeries
1D data vector in the form of a list/ numpy array or timeseries.
- i_roll: int
The number of indices to roll the array.
- Returns:
- rolled_waveform: 1d array or(pycbc TimeSeries object
The rolled wavefrom.
- waveformtools.waveformtools.save_obj(obj, name, obj_dir='./', protocol=5)[source]
A function to save python objects to disk using pickle.
- Parameters:
- objobject
The python object to be saved.
- namestring
The filename.
- obj_dirstring
The path to directory to be saved in. Defaults to PWD
- protocolint
The protocol to be used to save. Default is binary.
- Returns:
- Nothing (other than saving the data to the disk).
Notes
Protocols:
0 : Text 5 : Binary
See the man page of pickle for more details.
- waveformtools.waveformtools.shiftmatched(hdat, ind, delta_t=None, is_ts=False)[source]
Timeshift an array. IMP: After timeshifting, the original length of the array is retained by clipping last(first) when ind > 0(ind <0) ind number of points!!. Make sure the input array already has number of zeros z > ind (z<ind) initially at the end.
- Parameters:
- hdat: 1d array or a pycbc TimeSeries object
The input waveform to be shifted in time.
- ind: int
The numper of places to shift the input waveform.
- delta_t: int
the grid spacing in time.
- Returns:
- shifted_wf: a pycbc TimeSeries object
The waveform array of same length timeshifted by ind units by prepending zeros.
- waveformtools.waveformtools.shorten(tsdata, start, end, delta_t=None)[source]
Shorten an array given the start and end points.
- Parameters:
- tsdata: 1d array or a pycbc TimeSeries object
The waveform data.
- start: int
The start index of the data.
- end: int
The end index of the data.
- delta_t: float, optional
The time stepping.
- Returns:
- short_ts: a pycbc TimeSeries object
The shortened data, clipped before start and after end.
- waveformtools.waveformtools.simplematch_wfs_old(waveforms, delta_t=None)[source]
Simple match the given waveforms. Does not clip the waveforms at either ends.
- Parameters:
- waveforms: list
A list of pairs [waveform A, waveform B] of waveforms.
- delta_t: float, optional
The time stepping.
- Returns:
- match: list
A list of dicts [{ Aligned waveforms} , {match score (float), shift (number)}] containing the match information for all the input waveform pairs.
Notes
The time stepping delta_t is the same for each pair of waveforms in the list.
- waveformtools.waveformtools.smoothen(func_x, win, order, xdata=None, to_plot=False)[source]
Use the Savitzky-Golay Filter to smoothen the data. Show the plots if plot=`yes`.
- Parameters:
- func_x: 1d array
The y axis.
- win: int
Window for smoothening. Must be odd.
- order: int
The order of the polynomial used for interpolation.
- x: 1d array, optional.
The 1D list or numpy array, to plot the smoothened function. Only required if to_plot=True.
- to_plot: bool
True or False. Whether or not to display the plot.
- Returns:
- ydata: 1d array
The Savgol filtered list.
- waveformtools.waveformtools.startend(data)[source]
Identify the start and endpoints of the data.
- Parameters:
- data: 1d array or a pycbc TimeSeries object
The input waveform.
- Returns:
- start_index, end_index: int (2)
The pair of indices denoting the start and end points of an array
Notes
The starting and ending index of the non-zero part of the data is the identification criterion. Requires the data to be exactly zero outside a certain domain.
- waveformtools.waveformtools.taper(data, delta_t=1, zeros=150)[source]
A method to taper and append additional zeros at either ends, using the taper function of the pycbc TimeSeries object.
- Parameters:
- data: 1d array or a pycbc TimeSeries
The waveform data as list or numpy array or pycbc timeseries.
- delta_t: float
The timestepping.
- zeros: int
The number of zeros to be added.
- Returns:
- tapered_data: 1d array or a pycbc TimeSeries object
The waveform data tapered and zero padded. Returns a pycbc TimeSeries object.
Notes
See taper_timeseries from pycbc.waveform.utils for more details.
- waveformtools.waveformtools.taper_tanh(waveform, time_axis=None, delta_t=None, duration=10, sides='both')[source]
Taper a waveform with a \(tanh\) function at either ends
- Parameters:
- waveform: 1d array
A 1d array of waveform data.
- delta_t: float, optional
The time stepping delta_t. Optional if time_axis is given.
- percent: int, optional
The percent of data to taper. This is equally distributed on either sides of the array. Defaults to 10.
- sides: str
A string indicating which sides to taper. beg tapers the beginning, end tapers the end
and both tapers both the ends.
- Returns
- ——-
- time_axis, waveform: 1d array
The timeaxis and the waveform data of the tapered waveform.
- waveformtools.waveformtools.taperlengtheq(data_a, data_b, delta_t=None)[source]
Taper and equalize the lengths of two arrays.
- Parameters:
- data_a: list
The input waveform A.
- data_b: list
The input waveform B.
- delta_t: float
The time steping.
- Returns:
- equalized_signals: list
The Tapered, length equalized waveforms data_a and data_b, and a flag denoting which waveform was changed, a or b.
- waveformtools.waveformtools.totalmass(mass_ratio, chirp_mass)[source]
Find total mass from mass ratio and chirpmass.
- Parameters:
- mass_ratio: float
The mass ratio of the system.
- mchirp: float
The chirp mass of the system.
- Returns:
- total_mass: float
The total mass of the system.
- waveformtools.waveformtools.unwrap_phase(phi0)[source]
Unwrap the phase by finding turning points in phi0. Finding turning points for unwrapping arctan2 function
- Parameters:
- phi01darray
The wrapped phase which takes values in the range (0, 2pi).
- Returns
- ——-
- phic1darray
The unwrapped phase.
- waveformtools.waveformtools.xtract_camp(tsdata_p, tsdata_x, to_plot=False)[source]
Given real and imaginary parts of a complex timeseries, extract the amplitude of the complex data vector : (tsdata_p + i * tsdata_x)
- Parameters:
- tsdata_p: 1d array / a pycbc TimeSeries
- tsdata_x: 1d array / a pycbc TimeSeries
The plus and cross polarized components of the waveforms.
- to_plot: bool
True or False. Whether to plot the data or not
- Returns:
- camp: 1d array
The 1d array of extracted amplitudes of the waveform.
- waveformtools.waveformtools.xtract_camp_phase(tsdata_1, tsdata_2)[source]
Wrapper for extracting the amplitude and the phase of the complex vector.
- Parameters:
- tsdata_1: 1d array or pycbc TimeSeries object
- tsdata_2: 1d array or pycbc TimeSeries object
The two input waveforms as timeseries vectors, the plus and cross polarized components.
- delta_t: float
The timestepping delta_t.
- Returns:
- amplitude: 1d array
A list containing complex amplitude (list) and phase (list).
- waveformtools.waveformtools.xtract_cphase(tsdata_p, tsdata_x, delta_t=None, to_plot=False)[source]
Given real and imaginary parts of a complex timeseries, extract the phase of the waveform :arctan_(Img(data)/Re(data))
- Parameters:
- tsdata_p: 1d array / a pycbc TimeSeries
- tsdata_x: 1d array / a pycbc TimeSeries
The plus and cross polarized components of the waveforms.
- delta_t: float, optional
The time step. Overrides the timestep from pycbc TS object if given.
- to_plot: bool, optional
True or False. Whether to plot the data or not
- Returns:
- phic: 1d array
The 1d array of the phase of the waveform.