Template Waveform Function does not produce anything after being run

[amy5678]

I have been provided a python code as part of an assignment in which I've to create a waveform template by inputting my own values for the given variables. However, no matter where in the cell I define these, the function does not produce anything after being run.

The code is given by:

def make_template(m1,m2,fs,T,inv_psd,d=None,phic=None,tc=None):
    """
    This function generates a GW signal from a compact binary coalescence
    You must supply it with:
    m1 - the mass of the primary object (m1>m2) (units if solar mass)
    m2 - the mass of the secondary object (units of solar mass)
    fs - the sampling frequency (units of Hz)
    inv_psd - the inverse of the power spectral density (units of Hz)
    d [optional]- the distance to the source (units of megaparsecs)
    phic [optional] - the phase at coalescence (units of radians)
    tc [optional] - the time of coalescence from start (ubits of seconds)
    """
    N = int(T*fs)           # the total number of time samples
    dt = 1 / fs             # the sampling time (sec)
    df = 1.0/T              # the frequency resolution
    f_low = 20.0            # lowest frequency of waveform (Hz)
    if d is None:           # set the distance to 1 Mpc if not specified
        d = 1.0

    # make waveform
    n = N + 1
    while n>N/2:         # keep trying until the waveform safely fits inside the observation
        hp, _ = lalsimulation.SimInspiralChooseTDWaveform(
                    m1 * lal.MSUN_SI, m2 * lal.MSUN_SI,   # masses we comvert to SI units
                    0, 0, 0, 0, 0, 0,                     # spins we set to zero
                    d*1e6*lal.PC_SI,                      # the distance we convert to SI units
                    0, 0,            # inclination and phase
                    0, 0, 0,         # longAscNodes, eccentricity, and meanPerAno
                    dt,              # the time resolution
                    f_low,f_low,     # the min and reference frequencies
                    lal.CreateDict(),           # an empty placeholder
                    lalsimulation.IMRPhenomD)   # the waveform approximant 
        n = hp.data.data.size     # get the length of the generated waveform
        f_low += 1    # increment the low frequency start to reduce the signal length
   
    # pad the data and place the template in the middle 
    st = np.zeros(N)             # make an empty vector
    st[:n] = hp.data.data        # fill in the start of it with the signal
    idx = np.argmax(np.abs(st))  # find the time index with the maximum amplitude
    st = np.roll(st,-idx+N//2)   # roll the data to put the max amplitude in the centre

    # whiten the template using the given inverse PSD
    st_w = whiten(st,inv_psd,fs)

    # highpass the data - cut out all noise and signal power below 20 Hz 
    st_wb = highpass(st_w,20,fs)

    # put the template back and then shift the template to the requested time and phase
    f = df*np.arange(N//2 + 1)         # define frequency vector
    phase_t, phase_p = 0.0, 0.0        # initialise phase corrections
    if tc is not None:                 # if we set the time of coalescence
        st_wb = np.roll(st_wb,-N//2)   # undo the roll to the centre
        phase_t = 2.0j*np.pi*f*tc      # define phase correction for time shift
    if phic is not None:               # if we set the phase at coalescence
        phase_p = 1.0j*phic            # define phase shift

    # Fourier transform and apply the phase correction
    st_wbf = np.fft.rfft(st_wb)*np.exp(-phase_t - phase_p) 
    st_wb = np.fft.irfft(st_wbf)    # inverse transform back to the time domain

    # make a time vector
    t = np.arange(int(fs*T))/fs

    # return template
    return t, st_wb

# function to whiten data
def whiten(x, inv_psd, fs):
    """
    This function whitened a timeseries by converting to the frequency domain
    and dividing through by the sqrt of the power spectral density.
    You must supply it with:
    x - input timeseries
    inv_psd - the inverse power spetrcal density (units of Hz)
    fs - the sampling frequency (units of Hz)  
    """

    Nt = x.size                     # the length of the timeseries
    dt = 1.0/fs                          # the sampling time
    freqs = np.fft.rfftfreq(Nt, dt)      # the frequency vectore for this timeseries

    # whitening: transform to freq domain, divide by sqrt(psd), then transform back, 
    # taking care to get normalization right.
    dwindow = tukey(Nt, alpha=1./16)             # define a tukey window
    hf = np.fft.rfft(x*dwindow)                  # apply Fourier transform to windowed data
    norm = 1./np.sqrt(1./(dt*2))                 # define normalisation
    white_hf = hf * np.sqrt(inv_psd) * norm      # devide through by sqrt(PSD) and normalised
    white_ht = np.fft.irfft(white_hf, n=Nt)      # inverse Fourier transform back to time domain
    
    # return whitened timeseries
    return white_ht

def highpass(x,flow,fs):
    """
    This function applies a 4th order high pass butterworth filter to a
    timeseries.
    You must supply it with:
    x - the input timeseries
    flow - the filter high pass frequency cut-off (units of Hz)
    fs - the sampling frequency (units if Hz)
    """

    # We need to suppress the high frequency noise (no signal!) with some bandpassing:
    bb, ab = butter(4, flow*2./fs, btype='highpass')    # make the filter
    normalization = np.sqrt((fs/2.0-flow)/(fs/2.0))     # compute normalisation
    return filtfilt(bb, ab, x) / normalization          # return after applying the filter and normalisation

I have also been told in my task script that an example of this would be to run the following:

t, template = gw.make_template(20,12,4096,16,inv_psd)

Is this to imply that the function is built-in?

I hope this makes sense, any help is appreciated as it is due tomorrow evening!

[jefsummers]

Do you call the functions? If so, show that code as well

[amy5678]

Do you call the functions? If so, show that code as well

I don't have any other codes available - basically the small, one line of code near the bottom of my question is supposed to be an example of what I have to do. But, when i run that bit of code myself, absolutely nothing happens, even after defining all the variables.

[jefsummers]

The function call returns t and template, but you don't do anything to display them.