SIRENA functions

A B C D E F G H I J K L M N
O P Q R S T U V W X Y Z

Search functions by name at Index.

int addFirstRow(ReconstructInitSIRENA *reconstruct_init, fitsfile **inLibObject, double samprate, int runF0orB0val, gsl_vector *E, gsl_vector *PHEIGHT, gsl_matrix *PULSE, gsl_matrix *PULSEB0, gsl_matrix *MF, gsl_matrix *MFB0, gsl_matrix *COVAR, gsl_matrix *WEIGHT, gsl_matrix *PULSEMaxLengthFixedFilter)

Located in file: tasksSIRENA.cpp

This function writes the first row of the library (without intermediate AB-related values, because it would be necessary to have at least two rows=energies in the library). It also writes the FIXFILTT and FIXFILTF HDUs with the optimal filters in the time and frequency domain (real and imag parts) with fixed legnth and the PRCLOFWM HDU with the precalculated values for optimal filtering and :option:`mode`=*WEIGHTM*.

Members/Variables

ReconstructInitSIRENA **reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values). In particular, this function uses mode and the noise_spectrum in order to run calculus_optimalFilter().

fitsfile **inLibObject

FITS object containing information of the library FITS file

double samprate

Sampling rate

int runF0orB0val

If FilterMethod = F0 \Rightarrow runF0orB0val = 1. If FilterMethod = B0 \Rightarrow runF0orB0val = 0

gsl_vector *E

First energy to be included in the library

gsl_vector *PHEIGHT

Pulse height associated to the first energy to be included in the library

gsl_matrix *PULSE

Pulse template associated to the first energy to be included in the library

gsl_matrix *PULSEB0

Pulse template without baseline associated to the first energy to be included in the library

gsl_matrix *MF

Matched filter associated to the first energy to be included in the library

gsl_matrix *MFB0

Matched filter (baseline subtracted) associated to the first energy to be included in the library

gsl_matrix *COVAR

Covariance matrix associated to the first energy to be included in the library

gsl_matrix *WEIGHT

Weight matrix associated to the first energy to be included in the library

gsl_matrix *PULSEMaxLengthFixedFilter

Pulse template whose length is largeFilter associated to the first energy to be included in the library

int align(double samprate, gsl_vector **vector1, gsl_vector **vector2)

Located in file: tasksSIRENA.cpp

Based on [GilPita+05]

This function aligns vector1 with vector2 (by delaying or moving forward vector2) assuming that vector1 and vector2 are shifted replicas of the same function.

From the discrete function x[n] (n=0,...,N-1,N) and according to the time shifting property of the Fourier transform:

& x[n]    <------> X[f]\\
& x[n-m]  <------> X[f] exp(-j2\cdot\pi\cdot m/N)

If \mathit{Shift} = m then \mathit{PhaseDueToTheShift}= 2\pi m/N and thus, m = \mathit{PhaseDueToTheShift}\cdot N/(2\pi)

  1. Declare variables
  2. FFT of vector1
  3. FFT of vector2
  4. (Phases of the FFT_vector1 and FFT_vector2) *size/(2\pi)
  5. Shift between the input vectors
  6. shiftdouble into shiftint (because we are working with samples)
  7. Move forward or delay vector1 depending on positive or negative shift

Members/Variables

double samprate

Sampling rate

gsl_vector **vector1

GSL vector with input vector

gsl_vector **vector2

GSL with input vector which is delayed or moved forward to be aligned with vector1

int calculateEnergy(gsl_vector *vector, int pulseGrade, gsl_vector *filter, gsl_vector_complex *filterFFT, int runEMethod, int indexEalpha, int indexEbeta, ReconstructInitSIRENA *reconstruct_init, int domain, double samprate, gsl_vector *Pab, gsl_matrix *PRCLWN, gsl_matrix *PRCLOFWM, double *calculatedEnergy, int numlags, double *tstartNewDev)

Located in file: tasksSIRENA.cpp

This function calculates the energy of a pulse (vector) depending on the EnergyMethod, OFNoise, and the FilterDomain selected from input parameters.

  1. OPTFILT and NSD (= I2R or I2RALL, I2RNOL or I2RFITTED): Optimal filter = Wiener filter (see Optimal Filtering by using the noise spectral density)

Once the filter template has been created (filter or filterFFT), pulse height analysis is performed by aligning the template with a pulse and multiplying each point in the template by the corresponding point in the pulse. The sum of these products is the energy.

In the practice, the alignment of the pulse relative to the trigger is not completely accurate, so a number of n lags could be used in order to find the peak value of the energy. The n peak values are fitted to a parabola to find the most accurate energy (LagsOrNot) and a corrected starting time.

  1. OPTFILT and WEIGHTM (= I2R or I2RALL, I2RNOL or I2RFITTED) (see Optimal Filtering by using the noise weight matrix from noise intervals)
  1. WEIGHT and WEIGHTN (see Covariance matrices)

Members/Variables

gsl_vector *vector

Pulse whose energy has to be determined

int pulseGrade

Grade of the input pulse (to decide whether a full or only a rough estimation of energy is required).

gsl_vector *filter

Optimal filter in time domain

gsl_vector_complex *filterFFT

Optimal filter in frequency domain

int runEMethod
int indexEalpha

Index of the energy lower than the energy of the pulse which is being analyzed

int indexEbeta

Index of the energy higher than the energy of the pulse which is being analyzed

ReconstructInitSIRENA **reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values)

int domain
double samprate

Sampling rate in Hz

gsl_vector *Pab

PAB column in the library

gsl_vector *PRCLWN

PRCLx column in the library

gsl_vector *PRCLOFWM

OFWx column in the library

double *calculatedEnergy

Calculated energy in eV.

int numlags

Number of lags (if option:EnergyMethod = OPTFILT or I2R or I2RALL or I2RNOL or I2RFITTED and OFNoise = NSD)

double tstartNewDev

Addional deviation of the starting time (if LagsOrNot = 1)

int calculateIntParams(ReconstructInitSIRENA *reconstruct_init, int indexa, int indexb, double samprate, int runF0orB0val, gsl_matrix *modelsaux, gsl_matrix *covarianceaux, gsl_matrix *weightaux, gsl_vector *energycolumn, gsl_matrix **Wabaux, gsl_matrix **TVaux, gsl_vector **tEcolumn, gsl_matrix **XMaux, gsl_matrix **YVaux, gsl_matrix **ZVaux, gsl_vector **rEcolumn, gsl_matrix **Pabaux, gsl_matrix **Dabaux, gsl_matrix **PrecalWMaux, gsl_matrix **optimalfiltersabFREQaux, gsl_matrix **optimalfiltersabTIMEaux, gsl_matrix *modelsMaxLengthFixedFilteraux, gsl_matrix **PabMaxLengthFixedFilteraux)

Located in file: tasksSIRENA.cpp

This function calculates some intermediate scalars, vectors and matrices (WAB, TV, tE, XM, YV, ZV, rE, PAB and DAB) for the interpolation and covariance methods. See Covariance matrices reconstruction method.

  • Declare variables and allocate GSL vectors and matrices
  • Calculate intermediate scalars, vectors and matrices
  • Free allocated GSL vectors and matrices

Members/Variables

ReconstructInitSIRENA **reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values). In particular, this function uses PulseLength

int indexa

Lower index of the library to calculate the intermediate params (\alpha)

int indexb

Higher index of the library to calculate the intermediate params (\beta)

double samprate

Sampling rate

int runF0orB0val

If FilterMethod = F0 \Rightarrow runF0orB0val = 1. If FilterMethod = B0 \Rightarrow runF0orB0val = 0

gsl_matrix *modelsaux

GSL input matrix with model template

gsl_matrix *covarianceaux

GSL input matrix with weight matrix

gsl_matrix *weightaux

GSL input matrix with weight matrix

gsl_vector *energycolumn

GSL input vector with list of energies

gsl_matrix **WAB

Input/output intermediate parameter

gsl_matrix **TVaux

Input/output intermediate parameter

gsl_vector **tEcolumn

Input/output intermediate parameter

gsl_matrix **XMaux

Input/output intermediate parameter

gsl_matrix **YVaux

Input/output intermediate parameter

gsl_matrix **ZVaux

Input/output intermediate parameter

gsl_vector **rEcolumn

Input/output intermediate parameter

gsl_matrix **Pabaux

Input/output intermediate parameter

gsl_matrix **Dabaux

Input/output intermediate parameter

gsl_matrix **precalWMaux

Input/output intermediate parameter

gsl_matrix **optimalfiltersabFREQaux

Input/output intermediate parameter

gsl_matrix **optimalfiltersabTIMEaux

Input/output intermediate parameter

gsl_matrix *modelsMaxLengthFixedFilteraux

Input/output intermediate parameter

gsl_matrix **PabMaxLengthFixedFilteraux

Input/output intermediate parameter

int calculateTemplate(ReconstructInitSIRENA *reconstruct_init, PulsesCollection *pulsesAll, PulsesCollection *pulsesInRecord, double samprate, gsl_vector **pulseaverage, double *pulseaverageHeight, gsl_matrix **covariance, gsl_matrix **weight, int inputPulseLength, gsl_vector **pulseaverageMaxLengthFixedFilter)

Located in file: tasksSIRENA.cpp

This function calculates the template (PULSE column in the library) of non piled-up pulses. Just in case in the detection process some piled-up pulses have not been distinguished as different pulses, a pulseheights histogram is built. It uses the pulseheights histogram (built by using the PHEIGHT column of the library), Tstart and quality to select the non piled-up pulses.

  1. Declare and initialize variables

  2. Before building the histogram, select the pulseheights of the pulses well separated from other pulses and with quality = 0

  3. Create the pulseheights histogram

  4. Calculate the pulseaverage only taking into account the valid pulses:

    • Check if the pulse is piled-up or not
    • Non piled-up pulses \Rightarrow Average them
  5. Calculate covariance and weight matrices

  6. Free allocated GSL vectors

Members/Variables

ReconstructInitSIRENA **reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values). In particular, this function uses PulseLength and EnergyMethod

PulsesCollection *pulsesAll

Collection of pulses found in the previous records

PulsesCollection *pulsesInRecord

Collection of pulses found in the current record

gsl_vector **pulseaverage

GSL vector with the pulseaverage (template) of the non piled-up pulses

double *pulseaverageHeight

Height value of the pulseaverage

gsl_matrix **covariance

GSL matrix with covariance matrix

gsl_matrix **weight

GSL matrix with weight matrix (inverse of covariance matrix)

gsl_vector **pulseaverageMaxLengthFixedFilter

GSL vector with the pulseaverage (template) whose length is largeFilter of the non piled-up pulses

int calculus_optimalFilter(int TorF, int intermediate, int mode, gsl_vector *matchedfiltergsl, long mf_size, double samprate, int runF0orB0val, gsl_vector *freqgsl, gsl_vector *csdgsl, gsl_vector **optimal_filtergsl, gsl_vector **of_f, gsl_vector **of_FFT, gsl_vector_complex **of_FFT_complex)

Located in file: tasksSIRENA.cpp

See description also at optimal filter chapter

This function calculates the optimal filter for a pulse whose matched filter (normalized template) is provided as input parameter, matchedfiltergsl. An optimal filter is just a matched filter that has been adjusted based on the noise spectrum of the system.

It is assumed that all pulses are scaled versions of a template. In the frequency domain (as noise can be frequency dependent), the raw data can be expressed as P(f)=E\cdot S(f)+N(f), where S(f) is the normalized model pulse shape in the frequency domain, N(f) is the power spectrum of the noise and E is the scalar amplitude for the photon energy.

The second assumption is that the noise is stationary, i.e., it does not vary with time. The amplitude of each pulse can then be estimated by minimizing (weighted least-squares sense) the difference between the noisy data and the model pulse shape, being the \chi^2 condition to be minimized:

\chi^2 = \int \frac{(P(f)-E \cdot S(f))^2}{\langle\lvert N(f)\lvert ^2\rangle} df

In the time domain, the amplitude is the best weighted (optimally filtered) sum of the values in the pulse

E = k \int p(t)\cdot of(t)

where of(t) is the time domain expression of optimal filter which in frequency domain

OF(f) = \frac{S^*(f)}{\langle\lvert N(f)\lvert ^2\rangle}

and k is the normalization factor to give E in units of energy

k = \int \frac{S(f)\cdot S^{*}(f)}{\langle\lvert N(f)\lvert ^2\rangle} df

Steps:

  • FFT calculus of the matched filter (normalized template)

    • Declare variables
    • Complex FFT values for positive and negative frequencies
    • FFT calculus
    • Generation of the frequencies (positive and negative)
    • Magnitude and argument for positive and negative frequencies
    • Free allocated GSL vectors
  • N(f)

  • To divide MatchedFilter(f)/N^2(f) \Rightarrow MatchedFilter(f) and N(f) must have the same number of points

    • if (mf_size < freqgsl->size)

      • if ((freqgsl->size)%mf_size == 0) \Rightarrow Decimate noise samples
      • else \Rightarrow It is necessary to work only with the positive frequencies so as not to handle the f=0 \Rightarrow N(f) interpolation (interpolatePOS())
    • else if (mf_size > freqgsl->size) \Rightarrow Error: Noise spectrum must have more samples than pulse spectrum

    • else if (mf_size == freqgsl->size) \Rightarrow It is not necessary to do anything

  • OptimalFilter = MatchedFilter'(f)/N^2(f)

  • Calculus of the normalization factor

  • Inverse FFT (to get the expression of the optimal filter in time domain)

    • Complex OptimalFilter(f) \Rightarrow Taking into account magnitude MatchedFilter(f)/N^2(f) and phase given by MatchedFilter(f)
  • Free allocated GSL vectors

Members/Variables

int TorF

If FilterDomain = T \Rightarrow TorF = 0; If FilterDomain = F \Rightarrow TorF = 1

int intermediate

If intermediate = 0 \Rightarrow Do not write an intermediate file; If intermediate = 1 \Rightarrow Write an intermediate file

int mode

If mode = 0 \Rightarrow CALIBRATION run (library creation); If mode = 1 \Rightarrow RECONSTRUCTION run (energy determination)

gsl_vector *matchedfiltergsl

Matched filter associated to the pulse (in general, from the interpolation between two matched filters of the library)

long mf_size

Matched filter size (samples)

double samprate

Sampling rate

int runF0orB0val

If FilterMethod = F0 \Rightarrow runF0orB0val = 1. If FilterMethod = B0 \Rightarrow runF0orB0val = 0.

gsl_vector *freqgsl

Frequency axis of the current noise spectral density (input)

gsl_vector *csdgsl

Current noise spectral density (input)

gsl_vector **optimal_filtergsl

Optimal filter in time domain (output)

gsl_vector **of_f

Frequency axis of the optimal filter spectrum (output)

gsl_vector **of_FFT

Optimal filter spectrum (absolute values) (output)

gsl_vector_complex **of_FFT_complex

Optimal filter spectrum (complex values) (output)

int convertI2R(ReconstructInitSIRENA **reconstruct_init, TesRecord **record, gsl_vector **invector)

Located in file: tasksSIRENA.cpp

This funcion converts the current space into a quasi-resistance space (see Quasi Resistance Space for I2R, I2RALL, I2RNOL and I2RFITTED modes). The input invector filled in with current values is filled in here with I2R, I2RALL, I2RNOL or I2RFITTED quasi-resistances at the output.

  • Read R0, I0_START, ADUCNV and IMIN from the RecordFile by using the pointer (*reconstruct_init)->record_file_fptr
  • Conversion according to EnergyMethod = I2R
  • Conversion according to EnergyMethod = I2RALL
    • Read also RL, LFILTER and TTR from the RecordFile
  • Conversion according to EnergyMethod = I2RNOL
    • I2RALL neglecting the circuit inductance
  • Conversion according to EnergyMethod = I2RFITTED
  • Because in runEnergy() the record (TesRecord type) is used \Rightarrow The I2R, I2RALL, I2RNOL or I2RFITTED transformed record has to be used

Members/Variables

ReconstructInitSIRENA **reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values)

TesRecord **record

Structure containing the input ADC record

gsl_vector **invector

Input current (ADC) vector & output resistance (I2R or I2RALL or I2RNOL or I2RFITTED) vector

int createDetectFile(ReconstructInitSIRENA *reconstruct_init, double samprate, fitsfile **dtcObject, int inputPulselength)

Located in file: tasksSIRENA.cpp

This function creates an intermediate FITS file with some useful info (during the development phase) if the intermediate input parameter is set to 1.

The intermediate FITS file will contain 2 HDUs:

  • PULSES HDU will contain some info about the found pulses: TSTART, I0 (the pulse itself), TEND, TAURISE, TAUFALL and QUALITY
  • TESTINFO HDU will contain columns FILDER (the low-pass filtered and differentiated records) and THRESHOLD

If file exists \Rightarrow Check clobber for overwritting. If it does not, then create it.

Members/Variables

ReconstructInitSIRENA **reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values). In particular, this function uses some of its values to call the intermediate file (detectFile), to write input parameters info in the Primary extension of the intermediate file and also clobber

double samprate

Sampling rate (calculated from input FITS keyword DELTAT)

fitsfile dtcObject

Object which contains information of the intermediate FITS file (used also by writeTestInfo() and writePulses()).

int inputPulseLength

PulseLength input parameter

int createHisto(gsl_vector *invector, int nbins, gsl_vector **xhistogsl, gsl_vector **yhistogsl)

Located in file: tasksSIRENA.cpp

This function builds the histogram of the input vector.

  • Histogram x-axis values are the different input vector values (pulseheights)
  • Histogram y-axis values are the the number of cases per unit of the variable on the horizontal axis
  1. Declare variables
  2. It will work with the positive elements of the input vector \Rightarrow invectoraux2
  3. Check if all the values of invector are the same \Rightarrow Histogram of only one bin
  4. Obtain invector_max and invector_min
  5. Obtain binSize
  6. Create histogram axis
  7. Free allocated GSL vectors

Members/Variables

gsl_vector *invector

GSL input vector

int nbins

Number of bins to build the histogram

gsl_vector **xhistogsl

GSL vector with output histogram x-axis

gsl_vector **yhistogsl

GSL vector with output histogram y-axis

int createLibrary(ReconstructInitSIRENA *reconstruct_init, bool *appendToLibrary, fitsfile **inLibObject, int inputPulseLength)

Located in file: tasksSIRENA.cpp

This function creates the pulse templates library FITS file, if it does not exist yet. Otherwise, it opens it (to add a new row).

  1. If it exists \Rightarrow Open it and set appendToLibrary = true

  2. If it does not exist \Rightarrow Create it and set appendToLibrary = false

    • Write keyword EVENTCNT = 1 in the LIBRARY extension
    • Write info keywords about the running parameters in the Primary extension

Members/Variables

ReconstructInitSIRENA **reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values). In particular, this function uses parameters to call the library (LibraryFile) and to write input parameters info in the Primary extension of the library file

bool appendToLibrary

Used by the function writeLibrary()

fitsfile **inLibObject

Object which contains information of the library FITS file (used also by writeLibrary())

int inputPulseLength

PulseLength input parameter

int differentiate(gsl_vector **invector, int szVct)

Located in file: pulseprocess.cpp

This function applies the derivative method x_i-x_{i-1} to the input vector.

The derivative method provides more sensitivity to handle with piled-up pulses. Moreover, little variations of the baseline will not affect.

Members/Variables

gsl_vector **invector

Input/Ouput GSL vector (non-differentiate input vector/differentiate input vector)

int szVct

Size of invector

int eigenVV(gsl_matrix *matrixin, gsl_matrix **eigenvectors, gsl_vector **eigenvalues)

Located in file: tasksSIRENA.cpp

This funcion provides the principal eigenvectors and eigenvalues of the input matrix (at the moment, the first two eigenvalues and eigenvectors). The eigenvalues and eigenvectors are sorted in descending order and only the principal components are provided.

  • Calculate the eigenvectors and the eigenvalues
  • Sort the eigenvectors and the eigenvalues in descending order
  • Choose the main eigenvectors and eigenvalues (the principal components analysis). At the moment, the first two eigenvectors and eigenvalues

Members/Variables

gsl_matrix *matrixin

Input GSL matrix

gsl_matrix **eigenvectors

Subset of eigenvectors of ‘matrixin’ chosen by PCA (the first two ones)

gsl_vector **eigenvalues

Subset of eigenvalues of ‘matrixin’ chosen by PCA (the first two ones)

void exit_error(const char *const func, string msg, int status)

Located in file: genutils.cpp

This function prints out error messages and exits program.

Members/Variables

const char *const func

Function name whose error is printed

string msg

Error message to be printed

int status

Status

int FFT(gsl_vector *invector, gsl_vector_complex *outvector, double STD)

Located in file: genutils.cpp

This function calculates the FFT of the elements of a vector.

GSL library (overview of FFTs):

For physical applications it is important to remember that the index appearing in the DFT does not correspond directly to a physical frequency. If the time-step of the DFT is \Delta then the frequency domain includes both positive and negative frequencies, ranging from -1/(2\Delta) through 0 to +1/(2\Delta). The positive frequencies are stored from the beginning of the array up to the middle, and the negative frequencies are stored backwards from the end of the array.

Here is a table which shows the layout of the array data, and the correspondence between the time domain data z, and the frequency domain data x.

index z x = FFT(z)
0 z(t = 0) x(f = 0)
1 z(t = 1) x(f = 1/(n\Delta))
2 z(t = 2) x(f = 2/(n\Delta))
[…] [……..] [………………]
n/2 z(t = n/2) x(f = +1/(2\Delta),-1/(2\Delta))
[…] [……..] [………………]
n-3 z(t = n-3) x(f = -3/(n\Delta))
n-2 z(t = n-2) x(f = -2/(n\Delta))
n-1 z(t = n-1) x(f = -1/(n\Delta))

The frequency axis will be built as f = i/STD = i/(size/samprate) with i varying from 0 to size/2-1 (n=size and \Delta=1/samprate sec/sample).

Members/Variables

gsl_vector *invector

Input GSL vector

gsl_vector_complex *outvector

Output GSL complex vector with the FFT of invector

double STD

SelectedTimeDuration = (Size of invector)/samprate

int FFTinverse(gsl_vector_complex *invector, gsl_vector *outvector, double STD)

Located in file: genutils.cpp

This function calculates the inverse FFT of the elements of a vector.

Members/Variables

gsl_vector_complex *invector

Input GSL complex vector

gsl_vector *outvector

Output GSL vector with the inverse FFT of invector

double STD

SelectedTimeDuration = (Size of invector)/samprate

int filderLibrary(ReconstructInitSIRENA **reconstruct_init, double samprate)

Located in file: tasksSIRENA.cpp

This function calculates the (low-pass filter and) derivative of the models (pulse_templates) in the library (only necessary if first record), and it stores the pulse_templates_filder and the maxDERs and samp1DERs (see below) in the reconstruct_init structure.

The maximum of the (low-pass filtered and) differentiated pulse has to be compared to the maxDERs to select the appropriate model. Or, the 1st sample out of the differentiated pulse has to be compared to the samp1DERs to select the appropriate model.

  1. Check if it is the first record
  2. (Low-pass filter and) differentiate the models (pulse_templates) of the library
  3. Store the (low-pass filtered) derivatives in pulse_templates_filder
  4. Calculate the maximum of the (low-pass filtered and) differentiated models (maxDERs)
  5. Locate the 1st sample of the (low-pass filtered and) differentiated models (samp1DERs)

Members/Variables

ReconstructInitSIRENA **reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values). In particular, this function uses some parameters to filter (scaleFactor) and other parameters to handle the pulse templates and their derivatives (ntemplates, pulse_templates, pulse_templates_filder and maxDERs)

double samprate

Sampling rate

bool fileExists(const std::string &name)

Located in file: genutils.cpp

This function checks for file existence returning a boolean value.

Members/Variables

const std::string &name

File name

int filterByWavelets(ReconstructInitSIRENA *reconstruct_init, gsl_vector **invector, int length, int *onlyOnce)

Located in file: tasksSIRENA.cpp

This function filters the input/output signal invector, reducing the noise level.

Steps:

  • It is only going to work with n elements of :cpp:member:’invector’
  • Discrete Wavelet Treansform
  • Sorting coefficients
  • Hard thresholding: n-nc coefficeints are deleted (those with low energy)
  • Inverse DWT

Members/Variables

ReconstructInitSIRENA **reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values).

gsl_vector **invector

Input/output signal

int length

Length of the wavelet transform

int *onlyOnce

In order to control the times to be executed

int findMeanSigma(gsl_vector *invector, double *mean, double *sigma)

Located in file: pulseprocess.cpp

This function calculates the mean and the standard deviation of the input vector.

Members/Variables

gsl_vector *invector

Input GSL vector

double *mean

Mean of the elements of invector

double *sigma

Standard deviation of the elements of invector

int findPulsesCAL(gsl_vector *vectorin, gsl_vector *vectorinDER, gsl_vector **tstart, gsl_vector **quality, gsl_vector **pulseheight, gsl_vector **maxDERgsl, int *nPulses, double *threshold, double scalefactor, double samplingRate, int samplesup, double nsgms, double lb, double lrs, ReconstructInitSIRENA *reconstruct_init, double stopcriteriamkc, double kappamkc)

Located in file: pulseprocess.cpp

This function is going to find the pulses in a record (in the CALibration mode) by using the function findTstartCAL().

Steps:

  • Declare variables
  • Establish the threshold (call medianKappaClipping())
  • Find pulses (call findTstartCAL())
  • If at least a pulse is found
    • Get pulseheight of each found pulse (in order to be used to build the pulse templates library)
  • Free allocated GSL vectors

Members/Variables

gsl_vector *vectorin

Not filtered record

gsl_vector *vectorinDER

Derivative of the (low-pass filtered) vectorin

gsl_vector **tstart

Starting time of the found pulses into the record (in samples)

gsl_vector **quality

Quality of the found pulses into the record

gsl_vector **pulseheight

Pulse height of the found pulses into the record

gsl_vector **maxDERgsl

Maximum of the derivative of the found (low-pass filtered) pulses into the record

int *nPulses

Number of found pulses

double *threshold

Threshold used to find the pulses (output parameter because it is necessary out of the function)

double scalefactor

Scale factor to calculate the LPF box-car length (scaleFactor)

double samplingRate

Sampling rate

int samplesup

Number of consecutive samples over the threshold to locate a pulse (samplesUp)

double nsgms

Number of Sigmas to establish the threshold (nSgms)

double lb

Vector containing the baseline averaging length used for each pulse

double lrs

Running sum length (LrsT in samples)

ReconstructInitSIRENA *reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values). In particular, this function uses maxPulsesPerRecord (EventListSize) and PulseLength (and findTstartCAL() uses tstartPulse1, tstartPulse2 and tstartPulse3)

double stopcriteriamkc

Used in medianKappaClipping() (%)

double kappamkc

Used in medianKappaClipping()

int FindSecondaries(int maxPulsesPerRecord, gsl_vector *adjustedDerivative, double adaptativethreshold, ReconstructInitSIRENA *reconstruct_init, int tstartFirstEvent, int *numberPulses, gsl_vector **tstartgsl, gsl_vector **flagTruncated, gsl_vector **maxDERgsl, gsl_vector **lagsgsl)

Located in file: pulseprocess.cpp

This function runs after InitialTriggering() to find all the events (except the first one) in the first derivative of the (low-pass filtered) record by using the Adjusted Derivative detection method.

Steps:

  • Declare variables

  • Establishing the criteria of the slope of the derivative depending on the sampling rate

  • It is necessary to find the tstarts…

    It looks for an event and if a pulse is found, it looks for another event

    • It looks for an event since the beginning (or the previous event) to the end of the record. The first condition to detect an event is that the adjustedDerivative was over the threshold

      • Select the model of the found pulse from the libary by using the 1st sample of the derivative (samp1DER)

      • Dot product between the detected pulse and the pulse template in 3 different lags

        • If maximum of the dot product found \Rightarrow Stop calculating dot products in more lags
        • If maximum of the dot product not found \Rightarrow Calculate dot products in more lags (number of lags is limited to 5)
      • If maximum of the dot product not found \Rightarrow tstart is the first sample crossing above the threshold (without jitter)

        • Average of the first 4 samples of the derivative
        • Find model in order to subtract
      • If maximum of the dot product found \Rightarrow Parabola analytically defined \Rightarrow Locate the maximum \Rightarrow New tstart (with jitter)

        • Iterative process in order to extract the best template from the library:
          • samp1DER correction
          • Find the model from the libary by using the corrected samp1DER
          • Dot product in 3 lags
          • Locate the maximum of the parabola
        • samp1DER correction
        • Find model in order to subtract
        • Template correction
        • Average of the first 4 samples of the derivative
      • The second condition to detect an event is meeting the criteria of the slope of the derivative

    • Subtract the model from the adjusted derivative

      • Select the model of the found event from the libary by using the first sample of the derivative
      • Subtract
  • … Or to use the tstart provided as input parameters

    • Obtain the maxDERs of the events whose tstarts have been provided (by using the maximum of the derivative to find the model)
  • Free allocated GSL vectors

Members/Variables

int maxPulsesPerRecord

Expected maximum number of events per record in order to not allocate the GSL variables with the size of the input vector (EventListSize)

gsl_vector *adjustedDerivative

First derivative of the (low-pass filtered) record

double adaptativethreshold

Threshold

double samprate

Sampling rate

ReconstructInitSIRENA *reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values). In particular, this function uses some parameters (PulseLength, tstartPulse1, tstartPulse2, tstartPulse3, EnergyMethod) and the templates

int tstartFirstEvent

Tstart of the first event of the record (in samples) found by InitialTriggering()

int *numberPulses

Number of found events

gsl_vector **tstartgsl

Starting time of the found events (in samples)

gsl_vector **flagTruncated

Flag indicating if the event is truncated (inside this function only initial truncated pulses are classified)

gsl_vector **maxDERgsl

Maximum of the derivative of the event

gsl_vector **samp1DERgsl

Average of the first 4 samples of the derivative of the event

gsl_vector **lagsgsl

Number of necessary lags to establish the tstart (currently limited to 5)

int FindSecondariesSTC(int maxPulsesPerRecord, gsl_vector *adjustedDerivative, double adaptativethreshold, ReconstructInitSIRENA *reconstruct_init, int tstartFirstEvent, int *numberPulses, gsl_vector **tstartgsl, gsl_vector **flagTruncated, gsl_vector **maxDERgsl, gsl_vector **lagsgsl)

Located in file: pulseprocess.cpp

This function runs after InitialTriggering() to find all the events (except the first one) in the first derivative of the (low-pass filtered) record by using the Single Threshold Crossing method.

Steps:

  • Declare variables

  • It is necessary to find the tstarts…

    It looks for an event and if a pulse is found, it looks for another event

    • It looks for an event since the beginning (or the previous event) to the end of the record. The condition to detect an event is that the adjustedDerivative was over the threshold at least nSamplesUp consecutive samples
  • … Or to use the tstart provided as input parameters

    • Obtain the maxDERs of the events whose tstarts have been provided
  • Free allocated GSL vectors

Members/Variables

int maxPulsesPerRecord

Expected maximum number of events per record in order to not allocate the GSL variables with the size of the input vector (EventListSize)

gsl_vector *adjustedDerivative

First derivative of the (low-pass filtered) record

double adaptativethreshold

Threshold

double samprate

Sampling rate

ReconstructInitSIRENA *reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values). In particular, this function uses some parameters (PulseLength, tstartPulse1, tstartPulse2, tstartPulse3, EnergyMethod) and the templates

int tstartFirstEvent

Tstart of the first event of the record (in samples) found by InitialTriggering()

int *numberPulses

Number of found events

gsl_vector **tstartgsl

Starting time of the found events (in samples)

gsl_vector **flagTruncated

Flag indicating if the event is truncated (inside this function only initial truncated pulses are classified)

gsl_vector **maxDERgsl

Maximum of the derivative of the event

gsl_vector **samp1DERgsl

Average of the first 4 samples of the derivative of the event

int findTstartCAL(int maxPulsesPerRecord, gsl_vector *der, double adaptativethreshold, int nSamplesUp, ReconstructInitSIRENA *reconstruct_init, int *numberPulses, gsl_vector **tstartgsl, gsl_vector **flagTruncated, gsl_vector **maxDERgsl)

Located in file: pulseprocess.cpp

This function scans all the values of the derivative of the (low-pass filtered) record until it finds nSamplesUp consecutive values (due to the noise more than 1 value is required) over the threshold. To look for more pulses, it finds nSamplesUp consecutive values (due to the noise) under the threshold and then, it starts to scan again.

Steps:

  • Declare variables

  • Allocate GSL vectors

  • It is possible to find the tstarts…

    • Obtain tstart of each pulse in the derivative:
      • If der_i>threshold and foundPulse=false, it looks for nSamplesUp consecutive samples over the threshold
        • If not, it looks again for a pulse crossing over the threshold
        • If yes, a pulse is found (truncated if it is at the beginning)
      • If der_i>threshold and foundPulse=true, it looks for a sample under the threshold
        • If not, it looks again for a sample under the threshold
        • If yes, it looks for nSamplesUp consecutive samples under the threshold and again it starts to look for a pulse
  • … Or to use the tstart provided as input parameters

    Obtain the maxDERs of the pulses whose tstarts have been provided

Members/Variables

int maxPulsesPerRecord

Expected maximum number of pulses per record in order to not allocate the GSL variables with the size of the input vector (EventListSize)

gsl_vector *der

First derivative of the (low-pass filtered) record

double adaptativethreshold

Threshold

int nSamplesUp

Number of consecutive samples over the threshold to ‘find’ a pulse (samplesUp)

ReconstructInitSIRENA *reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values). In particular, this function uses PulseLength, tstartPulse1, tstartPulse2 and tstartPulse3

int *numberPulses

Number of found pulses

gsl_vector **tstartgsl

Pulses tstart (in samples)

gsl_vector **flagTruncated

Flag indicating if the pulse is truncated

gsl_vector **maxDERgsl

Maximum of the first derivative of the (low-pass filtered) record inside each found pulse

int find_Esboundary(double maxDER, gsl_vector *maxDERs, ReconstructInitSIRENA *reconstruct_init, int *indexEalpha, int *indexEbeta, double *Ealpha, double *Ebeta)

Located in file: tasksSIRENA.cpp.

This function provides the indexes of the two energies which straddle the pulse energy, by comparing the maximum value of the pulse derivative (maxDER) to the list of maximums in the library (maxDERs).

It finds the two embracing maxDERs in the calibration library

Members/Variables

double maxDER

Max value of the derivative of the (filtered) pulse whose embracing energies are being sought

gsl_vector *maxDERs

GSL vector with the maximum values of the derivatives of the templates in the library to be compared with the pulse being analysed

ReconstructInitSIRENA *reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values). In particular, this function uses the info in the library about the energies

int *indexEalpha

Index of the energy lower than the energy of the pulse which is being analyzed

int *indexEbeta

Index of the energy higher than the energy of the pulse which is being analyzed

double *Ealpha

Energy (in eV) which straddle the maxDER in the lower limit

double *Ebeta

Energy (in eV) which straddle the maxDER in the higher limit

int find_matchedfilter(int runF0orB0val, double maxDER, gsl_vector *maxDERs, ReconstructInitSIRENA *reconstruct_init, gsl_vector **matchedfilterFound, double *Ealpha, double *Ebeta)

Located in file: tasksSIRENA.cpp

This function chooses the proper matched filter (normalized template) of the calibration library (MF or MFB0) by comparing the maximum of their derivatives (maxDER versus maxDERs).

It finds the two embracing maxDERs in the calibration library and interpolates (interpolate_model())

  • If maxDER is lower than the lowest maxDERs in the library \Rightarrow The matched filter with the lowest maxDERs (first row) in the library is chosen
  • If maxDER is higher than the highest maxDERs in the library \Rightarrow The matched filter with the highest maxDERs (last row) in the library is chosen

Members/Variables

int runF0orB0val

If FilterMethod = F0 \Rightarrow runF0orB0val = 1. If FilterMethod = B0 \Rightarrow runF0orB0val = 0

double maxDER

Max value of the derivative of the (filtered) pulse whose matched filter is being sought

gsl_vector *maxDERs

GSL vector with the maximum values of the derivatives of the templates in the library to be compared with the pulse being analysed

ReconstructInitSIRENA *reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values). In particular, this function uses the info in the library (matched_filters)

gsl_vector **matchedfilterFound

GSL vector with the matched filter selected

double *Ealpha

Energy (in eV) which straddle the maxDER in the lower limit

double *Ebeta

Energy (in eV) which straddle the maxDER in the higher limit

int find_matchedfilterDAB(double maxDER, gsl_vector *maxDERs, ReconstructInitSIRENA *reconstruct_init, gsl_vector **matchedfilterFound, gsl_vector **PabFound, double *Ealpha, double *Ebeta)

Located in file: tasksSIRENA.cpp

This function selects the proper matched filter (normalized template) from the calibration library (from column DAB) by comparing the maximum value of the pulse derivative (maxDER) to the list of maximums in the library (maxDERs) for the DAB interpolation method (see optimal filter chapter). It also selects the proper row from the column PAB.

It finds the two embracing maxDERs in the calibration library

  • If maxDER is lower than the lowest maxDERs in the library \Rightarrow The data with the lowest maxDERs (first row) in the library are chosen
  • If maxDER is higher than the highest maxDERs in the library \Rightarrow The data of the penultimate row in the library are chosen

Members/Variables

int runF0orB0val

If FilterMethod = F0 \Rightarrow runF0orB0val = 1. If FilterMethod = B0 \Rightarrow runF0orB0val = 0

double maxDER

Max value of the derivative of the (filtered) pulse whose matched filter is being sought

gsl_vector *maxDERs

GSL vector with the maximum values of the derivatives of the templates in the library to be compared with the pulse being analysed

ReconstructInitSIRENA *reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values). In particular, this function uses the info in the library (matched_filters)

gsl_vector **matchedfilterFound

GSL vector with the matched filter selected

gsl_vector **PabFound

PAB column from the library

double *Ealpha

Energy (in eV) which straddle the maxDER in the lower limit

double *Ebeta

Energy (in eV) which straddle the maxDER in the higher limit

int find_model_energies(double energy, ReconstructInitSIRENA *reconstruct_init, gsl_vector **modelFound)

Located in file: pulseprocess.cpp

This function uses energy in order to choose the proper pulse template (pulse_templates_B0) of the calibration library.

In general, it finds the two energies wich straddle energy in the calibration library and interpolates (interpolate_model()):

  • If energy is lower than the lowest energy in the library \Rightarrow The model with the lowest energy in the library is chosen
  • If energy is higher than the highest energy in the library \Rightarrow The model with the highest energy in the library is chosen

Members/Variables

double energy

Energy of the pulse whose pulse template is being sought

ReconstructInitSIRENA *reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values). In particular, this function uses the energies of the models (energies) and their templates (pulse_templates), the number of templates in the library (ntemplates), the template duration (template_duration) and the pulse_templates_B0.

gsl_vector **modelFound

Found template of the pulse whose energy is energy

int find_model_maxDERs(double maxDER, ReconstructInitSIRENA *reconstruct_init, gsl_vector **modelFound)

Located in file: pulseprocess.cpp

This function uses the maximum of the derivative of the (filtered) pulse (maxDER) in order to choose the proper pulse template (pulse_templates_filder) of the calibration library.

In general, it finds the two maxDER which straddle maxDER in the calibration library and interpolates (interpolate_model()):

  • If maxDER is lower than the lowest maxDERs in the library \Rightarrow The model with the lowest maxDERs in the library is chosen
  • If maxDER is higher than the highest maxDERs in the library \Rightarrow The model with the highest maxDERs in the library is chosen

Members/Variables

double maxDER

Maximum of the derivative of the (filtered) pulse whose pulse template is being sought

ReconstructInitSIRENA *reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values). In particular, this function uses the number of templates in the library (ntemplates), the template duration (template_duration), the filtered and differentiated templates (pulse_templates_filder) and the maxDERs of the templates

gsl_vector **modelFound

Found template of the pulse whose maximum of the derivative of the filtered version is maxDER

int find_model_samp1DERs(double samp1DER, ReconstructInitSIRENA *reconstruct_init, gsl_vector **modelFound)

Located in file: pulseprocess.cpp

This function uses the 1st sample of the derivative of the filtered pulse (samp1DER) in order to choose the proper pulse template (pulse_templates_filder) of the calibration library.

It finds the two samp1DER closer in the calibration library and interpolates (interpolate_model())

  • If samp1DER is lower than the lowest samp1DER in the library \Rightarrow The model with the lowest samp1DER in the library is chosen
  • If samp1DER is higher than the highest samp1DER in the library \Rightarrow The model with the highest samp1DER in the library is chosen

Members/Variables

double samp1DER

1st sample of the derivative of the filtered pulse whose pulse template is being sought

ReconstructInitSIRENA *reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values). In particular, this function uses the 1st samples of the derivative of the models (samp1DERs) and their derived templates (pulse_templates_filder), the number of templates in the library (ntemplates) and the template duration (template_duration).

gsl_vector **modelFound

Found template of the pulse whose 1st sample of the derivative of the filtered pulse is samp1DER

int find_optimalfilter(double maxDER, gsl_vector *maxDERs, ReconstructInitSIRENA *reconstruct_init, gsl_vector **optimalfilterFound, double *Ealpha, double *Ebeta)

Located in file: tasksSIRENA.cpp

This function selects the proper optimal filter from the calibration library (columns Tx, ABTx, Fx or ABFx) by comparing the maximum value of the pulse derivative (maxDER) to the list of maximums in the library (maxDERs).

It finds the two embracing maxDERs in the calibration library and interpolates (interpolate_model())

  • If maxDER is lower than the lowest maxDERs in the library \Rightarrow The optimal filter with the lowest maxDERs (first row) in the library is chosen
  • If maxDER is higher than the highest maxDERs in the library \Rightarrow The optimal filter with the highest maxDERs (last row) in the library is chosen

Members/Variables

double maxDER

Max value of the derivative of the (filtered) pulse whose optimal filter is being sought

gsl_vector *maxDERs

GSL vector with the maximum values of the derivatives of the templates in the library to be compared with the pulse being analysed

ReconstructInitSIRENA *reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values). In particular, this function uses the info in the library (optimal_filters)

gsl_vector **optimalfilterFound

GSL vector with the optimal filter selected

double *Ealpha

Energy (in eV) which straddle the maxDER in the lower limit

double *Ebeta

Energy (in eV) which straddle the maxDER in the higher limit

int find_optimalfilterDAB(double maxDER, gsl_vector *maxDERs, ReconstructInitSIRENA *reconstruct_init, gsl_vector **optimalfilterFound, gsl_vector **PabFound, double *Ealpha, double *Ebeta)

Located in file: tasksSIRENA.cpp

This function selects the proper optimal filter from the calibration library (columns Tx, ABTx, Fx or ABFx) by comparing the maximum value of the pulse derivative (maxDER) to the list of maximums in the library (maxDERs). It also selects the proper row from the column PAB.

It finds the two embracing maxDERs in the calibration library

  • If maxDER is lower than the lowest maxDERs in the library \Rightarrow The data with the lowest maxDERs (first row) in the library are chosen
  • If maxDER is higher than the highest maxDERs in the library \Rightarrow The data of the penultimate row in the library are chosen

Members/Variables

double maxDER

Max value of the derivative of the (filtered) pulse whose optimal filter is being sought

gsl_vector *maxDERs

GSL vector with the maximum values of the derivatives of the templates in the library to be compared with the pulse being analysed

ReconstructInitSIRENA *reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values). In particular, this function uses the info in the library (optimal_filters)

gsl_vector **optimalfilterFound

GSL vector with the optimal filter selected

gsl_vector **PabFound

PAB column from the library

double *Ealpha

Energy (in eV) which straddle the maxDER in the lower limit

double *Ebeta

Energy (in eV) which straddle the maxDER in the higher limit

int find_prclofwm(double maxDER, gsl_vector *maxDERs, ReconstructInitSIRENA *reconstruct_init, gsl_vector **PRCLOFWMFound, double *Ealpha, double *Ebeta)

Located in file: tasksSIRENA.cpp

When EnergyMethod = OPTFILT and OFNoise = WEIGHTM this function selects the proper precalculated values (OFWx) from the calibration library by comparing the maximum value of the pulse derivative (maxDER) to the list of maximums in the library (maxDERs) for the OFLib =yes.

It finds the two embracing maxDERs in the calibration library

  • If maxDER is lower than the lowest maxDERs in the library \Rightarrow The data with the lowest maxDERs (first row) in the library are chosen
  • If maxDER is higher than the highest maxDERs in the library \Rightarrow The data of the penultimate row in the library are chosen

Members/Variables

double maxDER

Max value of the derivative of the (filtered) pulse whose optimal filter is being sought

gsl_vector *maxDERs

GSL vector with the maximum values of the derivatives of the templates in the library to be compared with the pulse being analysed

ReconstructInitSIRENA *reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values). In particular, this function uses the info in the library (PRCLOFWM)

gsl_vector **PRCLOFWMFound

GSL vector with some precalculated selected

double *Ealpha

Energy (in eV) which straddle the maxDER in the lower limit

double *Ebeta

Energy (in eV) which straddle the maxDER in the higher limit

int find_prclwn(double maxDER, gsl_vector *maxDERs, ReconstructInitSIRENA *reconstruct_init, gsl_vector **PRCLWNFound, gsl_vector **PabFound, double *Ealpha, double *Ebeta)

Located in file: tasksSIRENA.cpp

When EnergyMethod = WEIGHTN this function selects the proper precalculated values (PRCLx) from the calibration library by comparing the maximum value of the pulse derivative (maxDER) to the list of maximums in the library (maxDERs) for the OFLib =yes. It also selects the proper row from the column PAB.

It finds the two embracing maxDERs in the calibration library

  • If maxDER is lower than the lowest maxDERs in the library \Rightarrow The data with the lowest maxDERs (first row) in the library are chosen
  • If maxDER is higher than the highest maxDERs in the library \Rightarrow The data of the penultimate row in the library are chosen

Members/Variables

double maxDER

Max value of the derivative of the (filtered) pulse whose optimal filter is being sought

gsl_vector *maxDERs

GSL vector with the maximum values of the derivatives of the templates in the library to be compared with the pulse being analysed

ReconstructInitSIRENA *reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values). In particular, this function uses the info in the library (PRECALWN)

gsl_vector **PRCLWNFound

GSL vector with the precalculated values selected

gsl_vector **PabFound

PAB column from the library

double *Ealpha

Energy (in eV) which straddle the maxDER in the lower limit

double *Ebeta

Energy (in eV) which straddle the maxDER in the higher limit

extern_C_void freeOptimalFilterSIRENA(OptimalFilterSIRENA *OFilterColl)

Located in file: integraSIRENA.cpp

Destructor of OptimalFilterSIRENA structure.

Members/Variables

OptimalFilterSIRENA *OFilterColl

Instance of OptimalFilterSIRENA structure

extern_C_void freeReconstructInitSIRENA(ReconstructInitSIRENA *reconstruct_init)

Located in file: integraSIRENA.cpp

Destructor of ReconstructInitSIRENA structure.

Members/Variables

ReconstructInitSIRENA *reconstruct_init

Instance of ReconstructInitSIRENA structure

extern_C_void freePulsesCollection(PulsesCollection *PulsesColl)

Located in file: integraSIRENA.cpp

Destructor of PulsesCollection structure.

Members/Variables

PulsesCollection *PulsesColl

Instance of PulsesCollection structure

int fromGslMatrix(void **buffer, gsl_matrix **matrix, int type)

Located in file: inoututils.cpp

The function puts the values of the input GSL matrix into an output buffer.

Members/Variables

void **buffer

Output buffer

int fromGslVector(void **buffer, gsl_vector **array, int type)

Located in file: inoututils.cpp

The function puts the values of the input GSL vector into an output buffer.

Members/Variables

void **buffer

Output buffer

gsl_vector **array

Input GSL vector

. cpp:member:: int type

FITS type (TINT, TSHORT, TDOUBLE, etc.)
int getB(gsl_vector *vectorin, gsl_vector *tstart, int nPulses, gsl_vector **lb, int sizepulse, gsl_vector **B)

Located in file: pulseprocess.cpp

This function calculates the sum, B, of lb digitized data samples of a pulse-free interval immediately before each pulse. If the pulse-free interval before the current pulse is lower than lb, B is calculated with the available number of samples. If there is not a pulse-free interval before the pulse, it is looked for it after the current pulse. The number of samples of the pulse-free interval used to calculate B is stored in the lb vector.

Steps:

First of all, the auxiliary variable Baux is initialized to -999 and all the elements of the lb vector are equal to the LbT input parameter in samples. Then, the code is divided into 2 if statements:

  • When the current pulse is the first pulse into the record:
    • tstart \geq lb \Rightarrow Sum lb samples
    • 0<tstart<lb \Rightarrow Sum the available number of samples (although the available number of samples was lower than lb)
    • tstart=0 \Rightarrow If there is not a pulse-free interval before the pulse, it is looked for it after the current pulse
  • When the current pulse is not the first pulse into the record:
    • tstart_i-tend_{i-1} \geq lb \Rightarrow Sum lb samples
    • 0<tstart_i-tend{i-1}<lb \Rightarrow Sum the available number of samples (although the available number of samples was lower than lb)
    • If there is not a pulse-free interval before the pulse, it is looked for it after the current pulse

If Baux is still -999, a pulse-free interval can not be found to apply the running sum filter. This has to be taken into account, out of the function, to try to get a usable B.

Members/Variables

gsl_vector *vectorin

Input record

gsl_vector *tstart

Starting time of the pulses into the record

int nPulses

Number of pulses into the record

gsl_vector **lb

Vector containing the baseline averaging length used for each pulse

int sizepulse

Size of the pulse in samples

gsl_vector **B

In general, sum of the Lb digitized data samples (LbT input parameters in samples) of a pulse-free interval immediately before the current pulse

LibraryCollection *getLibraryCollection(const char *const filename, int mode, int hduPRECALWN, int hduPRCLOFWM, int largeFilter, char *filter_domain, int pulse_length, char *energy_method, char *ofnoise, char *filter_method, char oflib, char **ofinterp, double filtEev, int lagsornot, int *const status)

Located in file: integraSIRENA.cpp

This function creates and retrieves a LibraryCollection from a file.

  • Create LibraryCollection structure
  • Open FITS file in READONLY mode (move to the first HDU) and get number of templates (rows)
  • Allocate library structure
  • Get PULSE and MF column numbers (depending the different options)
  • Get template duration
  • Allocate library structure (cont.)
  • Get matched filter duration
  • Read different columns and populate the LibraryCollection structure
  • Added new code to handle the new HDUs FIXFILTF, FIXFILTT, PRECALWN and PRCLOFWM
  • Free allocated GSL vectors and matrices

Members/Variables

int mode

Calibration run (0) or energy reconstruction run (1), mode

char *energy_method

Energy calculation Method: OPTFILT, WEIGHT, WEIGHTN, I2R, I2RALL, I2RNOL or I2RFITTED, EnergyMethod

int hduPRECALWN

Add or not the PRECALWN HDU in the library file (1/0) (only for library creation, mode = 0), hduPRECALWN

int hduPRCLOFWM

Add or not the hduPRCLOFWM HDU in the library file (1/0) (only for library creation, mode = 0), hduPRCLOFWM

int largeFilter

Length of the longest fixed filters (only for library creation, mode = 0), largeFilter

char *filter_domain

Filtering Domain: Time (T) or Frequency (F), FilterDomain

int pulse_length

Pulse length, PulseLength

char *energy_method

Energy calculation Method: OPTFILT, WEIGHT, WEIGHTN, I2R, I2RALL, I2RNOL, I2RFITTED or PCA, EnergyMethod

char *ofnoise

For optimal filtering : NSD or WEIGHTM, OFNoise

char *filter_method

Filtering Method: F0 (deleting the zero frequency bin) or B0 (deleting the baseline), FilterMethod

char **ofinterp

Optimal Filter by using the Matched Filter or the DAB as matched filter (MF/DAB) It has been fixed in ‘tesreconstruction as DAB (but it would be possible to work with MF)

double filtEev

Energy of the filters of the library to be used to calculate energy (only for OPTFILT, I2R, I2RALL, I2RNOL and I2RFITTED), filtEev

int *const status

Input/output status

NoiseSpec *getNoiseSpec(const char *const filename, int mode, int hduPRCLOFWM, char *energy_method, char *ofnoise, char *filter_method, int *const status)

Located in file: integraSIRENA.cpp

This function creates and retrieves a NoiseSpec from a file.

  • Create NoiseSpec structure
  • Open FITS file, move to the NOISE, NOISEALL and WEIGHTMS HDUs and get interesting keywords
  • Allocate NoiseSpec structure
  • Get noise spectrum (CSD), and noise frequencies (FREQ) column numbers
  • Read column CSD and save it into the structure
  • Read column FREQ and save it into the structure
  • Read columns Wx with the noise weight matrix from noise intervals and save them into the structure
  • Return noise spectrum

Members/Variables

const char *const filename

File name with noise

int mode

Calibration run (0) or energy reconstruction run (1), mode

int hduPRCLOFWM

Add or not the hduPRCLOFWM HDU in the library file (1/0) (only for library creation, mode = 0), hduPRCLOFWM

char *energy_method

Energy calculation Method: OPTFILT, WEIGHT, WEIGHTN, I2R, I2RALL, I2RNOL, I2RFITTED or PCA, EnergyMethod

char *ofnoise

For optimal filtering: NSD or WEIGHTM, OFNoise

char *filter_method

Filtering Method: F0 (deleting the zero frequency bin) or B0 (deleting the baseline), FilterMethod

int *const status

Input/Output status

int getPulseHeight(gsl_vector *vectorin, double tstart, double tstartnext, int lastPulse, double lrs, double lb, double B, int sizepulse, double *pulseheight)

Located in file: pulseprocess.cpp

This function estimates the pulse height of a pulse by using a running sum filter. It extracts from the record, vectorin, the pulse whose pulse height is going to be estimated by using RS_filter().

Steps:

  • Declare variables
  • Extracting from the record the pulse whose pulse height is going to be estimated
  • Apply the running sum filter

Members/Variables

gsl_vector *vectorin

Not filtered record

double tstart

Starting time of the pulse whose pulse height is going to be estimated

double tstartnext

Starting time of the next pulse whose pulse height is going to be estimated

int lastPulse

It is 1 if the pulse is the last one into the record or the only one

double lrs

Running sum length (equal to the LrsT input parameter in samples)

double lb

Baseline averaging length used for the pulse whose pulse height is going to be estimated

double B

In general, sum of the Lb digitized data samples (LbT input parameters in samples) of a pulse-free interval immediately before the current pulse

int sizepulse

Size of the pulse in samples

double *pulseheight

Estimated pulse height of the pulse

void gsl_vector_complex_absIFCA(gsl_vector *cvnew, gsl_vector_complex *cv)

Located in file: genutils.cpp

This function calculates the magnitude of the complex elements of a vector (real part).

Members/Variables

gsl_vector_complex *cv

Input GSL complex vector

gsl_vector *cvnew

Output GSL vector with the absolute values of the elements of cv

void gsl_vector_complex_argIFCA(gsl_vector *varg, gsl_vector_complex *vin)

Located in file: genutils.cpp

This function calculates the arguments of the complex elements of a vector.

Members/Variables

gsl_vector_complex *vin

Input GSL complex vector

gsl_vector *varg

Output GSL vector with the arguments of the elements of vin

void gsl_vector_complex_scaleIFCA(gsl_vector_complex *cv, gsl_complex z)

Located in file: genutils.cpp

This function multiplies the complex elements of a vector by a complex number.

Members/Variables

gsl_vector_complex *cv

Input/Output (scaled) GSL complex vector

gsl_complex z

Input GSL complex number

void gsl_vector_sqrtIFCA(gsl_vector *cvnew, gsl_vector *cv)

Located in file: genutils.cpp

This function calculates the square root of the elements of a vector.

Members/Variables

gsl_vector *cv

Input GSL complex vector

gsl_vector *cvnew

Output GSL vector with the square root values of the elements of cv

int gsl_vector_Sumsubvector(gsl_vector *invector, long offset, long n, double *sum)

Located in file: genutils.cpp

This function returns the sum of some elements of the input vector.

The starting element of the sum is offset from the start of the input vector. It will sum up n elements.

offset can take values from 0 to invector->size

Members/Variables

gsl_vector *invector

Input GSL vector

long offset

It is the first element to be summed

long n

Number of elements in the sum

double *sum

Calculated output value (sum of the corresponding elements)

extern_C_void initializeReconstructionSIRENA(ReconstructInitSIRENA *reconstruct_init, char *const record_file, fitsfile *fptr, char *const library_file, char *const event_file, int pulse_length, double scaleFactor, double samplesUp, double samplesDown, double nSgms, int detectSP, int mode, char *detectionMode, double LrsT, double LbT, char *const noise_file, char *filter_domain, char *filter_method, char *energy_method, double filtEev, char *ofnoise, int lagsornot, int ofiter, char oflib, char *ofinterp, char *oflength_strategy, int oflength, double monoenergy, char hduPRECALWN, char hduPRCLOFWM, int largeFilter, int interm, char *const detectFile, char *const filterFile, char clobber, int maxPulsesPerRecord, double SaturationValue, int tstartPulse1, int tstartPulse2, int tstartPulse3, double energyPCA1, double energyPCA2, char *const XMLFile, int *const status)

Located in file: integraSIRENA.cpp

This function initializes the structure ReconstructInitSIRENA with the variables required for SIRENA reconstruction. The values are taken from the input parameters.

  • Load LibraryCollection structure if library file exists
  • Load NoiseSpec structure
  • Fill in reconstruct_init

Members/Variables

ReconstructInitSIRENA *reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values)

char *const record_file

Filename of input data file with records, RecordFile

fitsfile *fptr

FITS object with pointer to data file

char *const library_file

File name of calibration library, LibraryFile

char *const event_file

File name of output events (with reconstructed energy), TesEventFile

int pulse_length

Pulse length, PulseLength

double scaleFactor

Detection scale factor for initial filtering, scaleFactor

double samplesUp

Number of samples for threshold trespassing, samplesUp

double samplesDown

Number of samples below the threshold to look for other pulse, samplesDown

double nSgms

Number of standard deviations in the kappa-clipping process for threshold estimation, nSgms

int detectSP

Detect secondary pulses (1) or not (0), detectSP

int mode

Calibration run (0) or energy reconstruction run (1), mode

char *detectionMode

Adjusted Derivative (AD) or Single Threshold Crossing (STC), detectionMode

double LrsT

Running sum length for the RS raw energy estimation (seconds), LrsT

double LbT

Baseline averaging length for the RS raw energy estimation (seconds), LbT

char *const noise_file

Noise file, NoiseFile

char *filter_domain

Filtering Domain: Time (T) or Frequency (F), FilterDomain

char *filter_method

Filtering Method: F0 (deleting the zero frequency bin) or F0 (deleting the baseline), FilterMethod

char *energy_method

Energy calculation Method: OPTFILT, WEIGHT, WEIGHTN, I2R, I2RALL, I2RNOL, I2RFITTED or PCA, EnergyMethod

double filtEev

Energy of the filters of the library to be used to calculate energy (only for OPTFILT, I2R, I2RALL, I2RNOL and I2RFITTED), filtEev

char *ofnoise

For optimal filtering: NSD or WEIGHTM, OFNoise

int lagsornot

Lags (1) or no lags (0), LagsOrNot

int ofiter

Iterate (1) or not iterate (0), OFIter

char oflib

Work or not with a library with optimal filters (1/0), OFLib

char *ofinterp

Optimal Filter by using the Matched Filter or the DAB as matched filter (MF/DAB) It has been fixed in tesreconstruction as DAB

char *oflength_strategy

Optimal Filter length Strategy: FREE, BASE2, BYGRADE or FIXED, OFStrategy

int oflength

Optimal Filter length (taken into account if OFStrategy = FIXED), OFLength

double monoenergy

Monochromatic energy of input file in eV (only for library creation, mode = 0), monoenergy

int hduPRECALWN

Add or not the PRECALWN HDU in the library file (1/0) (only for library creation, mode = 0), hduPRECALWN

int hduPRCLOFWM

Add or not the hduPRCLOFWM HDU in the library file (1/0) (only for library creation, mode = 0), hduPRCLOFWM

int largeFilter

Length of the longest fixed filters (only for library creation, mode = 0), largeFilter

int interm

Write or not intermediate files (1/0), intermediate

char *const detectFile

Intermediate detections file (if intermediate = 1), detectFile

char *const filterFile

Intermediate filters file (if intermediate = 1), filterFile

char clobber

Overwrite or not output files if exist (1/0), clobber

int maxPulsesPerRecord

Default size of the event list, EventListSize

double SaturationValue

Saturation level of the ADC curves

int tstartPulse1

Tstart (samples) of the first pulse (different from 0 if the tstartPulsei input parameters are going to be used), tstartPulse1

int tstartPulse2

Tstart (samples) of the second pulse, tstartPulse2

int tstartPulse3

Tstart (samples) of the third pulse (if 0 \Rightarrow PAIRS, if not 0 \Rightarrow TRIOS), tstartPulse3

double energyPCA1

First energy (only for EnergyMethod = PCA)

double energyPCA2

Second energy (only for EnergyMethod = PCA)

char *const XMLFile

File name of the XML input file with instrument definition

int *const status

Input/Output status

int InitialTriggering(gsl_vector *derivative, double nSgms, double scalefactor, double samplingRate, double stopcriteriamkc, double kappamkc, bool *triggerCondition, int *tstart, int *flagTruncated, double *threshold, int tstartProvided)

Located in file: pulseprocess.cpp

This function finds the first pulse in the input vector, first derivative of the (low-pass filtered) record.

Steps:

  • Declare variables

  • Stablish the threshold

  • It is necessary to find the tstart of the first pulse…

    Obtain tstart of the first pulse in the derivative if derivative_i>threshold

  • … Or to use the tstart provided as input parameter

Members/Variables

gsl_vector *derivative

First derivative of the (low-pass filtered) record

double nSgms

Number of Sigmas to establish the threshold (nSgms)

double scalefactor

Scale factor to calculate the LPF box-car length (scaleFactor)

double samplingRate

Sampling rate

double stopcriteriamkc

Used in medianKappaClipping() (%)

double kappamkc

Used in medianKappaClipping()

bool *triggerCondition

True \Rightarrow The algorithm has found the first event

False \Rightarrow The algorithm has not found any event

int *tstart

First event tstart (in samples)

int *flagTruncated

Flag indicating if the event is truncated

double *threshold

Calculated threshold (output parameter because it is necessary out of the function)

int tstartProvided

Tstart of the first pulse provided as input parameter

int interactivePars(inparam *taskPars, int np, string task)

Located in file inoututils.cpp

This function reads input parameters interactively (provided by the user or taken as default values). Used in tool gennoisespec.

Members/Variables

inparam *taskPars

Instance of inparam structure storing input parameters

int np

Number of parameters

string task

Tool name

int interpolatePOS(gsl_vector *x_in, gsl_vector *y_in, long size, double step, gsl_vector **x_out, gsl_vector **y_out)

Located in file: tasksSIRENA.cpp

This function interpolates an input vector (x_in, y_in), creating an output vector (x_out, y_out) with the size and frequency step given. POS comes from the fact that the input spectrum only has positive frequencies (in order to not handle the f=0 bin).

  • Declare and initialize variables
  • GSL method applied for interpolatation
  • Generate the interpolated output vector
  • Free memory

Members/Variables

gsl_vector *x_in

GSL input vector with the abscissas of the vector which is going to be interpolated

gsl_vector *y_in

GSL input vector with the ordinates of the vector which is going to be interpolated

long size

Size of the interpolated output vector

double step

Frequency step of the interpolated output vector

gsl_vector **x_out

GSL output vector with the abscissas of the interpolated vector

gsl_vector **y_out

GSL output vector with the ordinates of the interpolated vector

int interpolate_model(gsl_vector **modelFound, double p_model, gsl_vector *modelIn1, double p_modelIn1, gsl_vector *modelIn2, double p_modelIn2)

Located in file: pulseprocess.cpp

This function interpolates the pulse model, p(t,E), between two models of the pulse models library, p(t,E_1) and p(t,E_2), being E_1<E<E_2.

According to the interpolation method:

p(t,E)={\frac{E_2-E}{E_2-E_1}}p(t,E_1)+{\frac{E-E_1}{E_2-E_1}}p(t,E_2)

Members/Variables

gsl_vector **modelFound

Found model of the pulse whose energy or maxDER is p_model

double p_model

Parameter (energy or maxDER) of the pulse whose model is being sought

gsl_vector *modelIn1

Model of the pulse whose parameter (energy or maxDER) is immediately lower than p_model in the library FITS file

double p_modelIn1

Parameter (energy or maxDER) immediately lower than p_model in the library FITS file

gsl_vector *modelIn2

Model of the pulse whose parameter (energy or maxDER) is immediately greater than p_model in the library FITS file

double p_modelIn2

Parameter (energy or maxDER) immediately greater than p_model in the library FITS file

int loadRecord(TesRecord *record, double *time_record, gsl_vector **adc_double)

Located in file: tasksSIRENA.cpp

This fucntion loads the structure record into the adc_double GSL vector.

It checks if the record has been filled out with 0’s \Rightarrow It only loads the first values (which are different from 0).

Members/Variables

TesRecord *record

Member of TesRecord structure that contains the input record

double time_record

Starting time of the record (output)

gsl_vector **adc_double

Storage of the record to be processed (input/output)

int lpf_boxcar(gsl_vector **invector, int szVct, int sampleRate)

Located in file: pulseprocess.cpp

This function implements a low pass filtering as a box-car function in time.

The box-car function is a temporal average window:

x_{i-1}=\sum_{0}^{n-1}\frac{I_i}{n}

x_i=\sum_{1}^{n}\frac{I_i}{n}

If the cut frequency of the filter is \mathit{f_c}, the box-car length (n) is

\frac{1}{f_c}samprate

Steps:

  • Declare variables
  • Define the LPF (frequency domain) and the box-car function (time domain)
  • It is going to work with a longer vector to not have fake results for the last boxLength windows
  • Apply the box-car window by shifting it along the (lengthened) input vector
  • Free allocated GSL vectors

The function returns:

  • 1: Function cannot run
  • 3: Cut-off frequency too high \Rightarrow Equivalent to not filter
  • 4: Cut-off frequency too low

Members/Variables

gsl_vector **invector

Input/Output GSL vector (non-filtered input vector/filtered input vector)

int szVct

Size of invector

int sampleRate

Sampling rate (samples/s)

int matrix2vector(gsl_matrix *matrixin, gsl_vector **vectorout)

Located in file: tasksSIRENA.cpp

This function converts an input [n \times n] square matrix into an output n^2 vector. It puts the first row of the matrix (n elements) in the first n elements of the vector (from 0 to n-1), the second row of the matrix in the elements from n to 2n-1 of the vector and so on.

Members/Variables

gsl_matrix *matrixin

GSL input square matrix [n \times n]

gsl_vector **vectorout

GSL output vector whose length is n^2

int medianKappaClipping(gsl_vector *invector, double kappa, double stopCriteria, double nSigmas, int boxLPF, double *threshold)

Located in file: pulseprocess.cpp

This function calculates a threshold in the first derivative of the record by using a Kappa-clipping method (replacing points beyond mean\pm kappa \cdot sigma with the median).

Mean and sigma are calculated and values of invector out of (mean+kappa \cdot sigma,mean-kappa \cdot sigma) are replaced with the median (it is trying to look for the baseline). And this process is iteratively repeated until there are no points beyond mean \pm kappa \cdot sigma. Finally, the threshold is calculated as mean+nSigmas \cdot sigma (‘+’ is used because if there are pulses in the input invector they are always positive).

Steps:

  • Declare variables
  • Calculate the median
  • Iterate until there are no points out of the maximum excursion ( kappa \cdot sigma)
  • Establish the threshold as mean+nSigmas*sigma

Members/Variables

gsl_vector *invector

First derivative of the (filtered) record

double kappa

Value to establish the range around of the mean

double stopCriteria

It is given in %

double nSigmas

Times sigma to calculate threshold as mean+nSigmas \cdot sigma

int boxLPF

Length of the low-pass filtering box-car

double *threshold

Calculated threshold

extern_C_OptimalFilterSIRENA *newOptimalFilterSIRENA(int *const status)

Located in file: integraSIRENA.cpp

Constructor. It returns a pointer to an empty OptimalFilterSIRENA data structure.

Members/Variables

int *const status

Input/output status

extern_C_PulsesCollection *newPulsesCollection(int *const status)

Located in file: integraSIRENA.cpp

Constructor. It returns a pointer to an empty PulsesCollection data structure.

Members/Variables

int *const status

Input/output status

extern_C_ReconstructInitSIRENA *newReconstructInitSIRENA(int *const status)

Located in file integraSIRENA.cpp

Constructor. It returns a pointer to an empty ReconstructInitSIRENA data structure.

Members/Variables

int *const status

Input/output status

int noDetect(gsl_vector *der, ReconstructInitSIRENA *reconstruct_init, int *numberPulses, gsl_vector **tstartgsl, gsl_vector **flagTruncated, gsl_vector **maxDERgsl, gsl_vector **samp1DERgsl)

Located in file: pulseprocess.cpp

This function runs if the starting time of the pulses are agiven as input parameters (tstartPulse1 != 0). It looks for the maximum of the derivative of the pulse and the average of the first 4 samples of the derivative of the pulse.

Members/Variables

gsl_vector *der

First derivative of the (low-pass filtered) record

ReconstructInitSIRENA *reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values). In particular, this function uses some parameters (tstartPulse1, tstartPulse2, tstartPulse3) and the templates

int *numberPulses

Number of events

gsl_vector **tstartgsl

Starting time of the events (in samples)

gsl_vector **flagTruncated

Flag indicating if the event is truncated (inside this function only initial truncated pulses are classified)

gsl_vector **maxDERgsl

Maximum of the derivative of the event

gsl_vector **samp1DERgsl

Average of the first 4 samples of the derivative of the event

int polyFit(gsl_vector *x_fit, gsl_vector *y_fit, double *a, double *b, double *c)

Located in file: genutils.cpp

This function makes a polynomial fitting ax^2+bx+c using the regression quadratic analysis. To measure how well model agrees with the data, the chi-square merit function is used, which in this case is

\chi^2 (a,b,c)= \sum_{i=1}^{N}\left(\frac{y_i-a-bx_i-c{x_i}^2}{\sigma_i}\right)^2

This equation is minimized to determine a, b and c. Then

\begin{array}{ccc}
S_{(x,x)}=\sum{{x_i}^2}-\frac{(\sum{x_i})^2}{N} & S_{(x^2,y)}=\sum{{x_i}^2y_i}-\frac{\sum{{x_i}^2}\cdot\sum{y_i}}{N}\\
S_{(x,y)}=\sum{x_iy_i}-\frac{\sum{x_i}\cdot\sum{y_i}}{N} & S_{(x^2,x^2)}=\sum{{x_i}^4}-\frac{\left(\sum{{x_i}^2}\right)^2}{N}\\
S_{(x,x^2)}=\sum{{x_i}^3}-\frac{\sum{x_i}\cdot\sum{{x_i}^2}}{N} & \\
\end{array}

a = \frac{S_{(x^2,y)}S_{(x,x)}-S_{(x,y)}S_{(x,x^2)}}{S_{(x,x)}S_{(x^2,x^2)} -{\vert S_{(x,x^2)} \vert}^2}

b = \frac{S_{(x,y)}S_{(x^2,x^2)}-S_{(x^2,y)}S_{(x,x^2)}}{S_{(x,x)}S_{(x^2,x^2)} -{\vert S_{(x,x^2)} \vert}^2}

c = \frac{\sum{y_i}}{N}-b\frac{\sum{x_i}}{N}-a\frac{\sum{{x_i}^2}}{N}

Members/Variables

gsl_vector *x_fit

Input GSL with x vector

gsl_vector *y_fit

Input GSL with y vector

double *a

Fit coefficient of the quadratic term

double *b

Fit coefficient of the linear term

double *c

Fit coefficient (independent term)

int polyFitLinear(gsl_vector *x_fit, gsl_vector *y_fit, double *a, double *b)

Located in file: genutils.cpp

This function makes a linear fitting ax+b using the regression linear analysis. To measure how well model agrees with the data, the chi-square merit function is used, which in this case is

\chi^2 (a,b)= \sum_{i=1}^{N}\left(\frac{y_i-a-bx_i}{\sigma_i}\right)^2

This equation is minimized to determine a and b. Then

a = \frac{N \sum{x_i y_i}- \sum{x_i}\sum{y_i}}{N \sum{x_i^2} - {\left(\sum{x_i}\right)}^2}

b = \frac{\sum{y_i}}{N}-a \frac{\sum{x_i}}{N}

Members/Variables

gsl_vector *x_fit

Input GSL with x vector

gsl_vector *y_fit

Input GSL with y vector

double *a

Fit coefficient of the linear term

double *b

Fit coefficient (independent term)

void print_error(const char *const func, string message, int status)

Located in file: genutils.cpp

This function prints out error messages.

Members/Variables

const char *const func

Function name whose error is printed

string msg

Error message to be printed

int status

Status

int procRecord(ReconstructInitSIRENA **reconstruct_init, double tstartRecord, double samprate, fitsfile *dtcObject, gsl_vector *record, PulsesCollection *foundPulses)

Located in file: tasksSIRENA.cpp

This function processes the input record:

  1. Declare and initialize variables
  2. Allocate GSL vectors
  3. (Low-pass filtering and) differentiation
  4. Find the events (pulses) in the record
  5. Calculate the end time of the found pulses and check if the pulses are saturated
  6. Obtain the approximate rise and fall times of each pulse (to be done)
  7. Load the found pulses data in the input/output foundPulses structure
  8. Write test info (if reconstruct_init->intermediate = 1)
  9. Write pulses info in intermediate output FITS file (if reconstruct_init->intermediate = 1)
  10. Free allocated GSL vectors

Members/Variables

ReconstructInitSIRENA **reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values). In particular, this function uses parameters to filter (scaleFactor), to find pulses (PulseLength, samplesUp, nSgms, LrsT, LbT, EventListSize) and to write info (detectFile) if intermediate = 1

double tstartRecord

Starting time of the record (in order to calculate absolute times)

double samprate

Sampling rate (in order to low-pass filter)

fitsfile *dtcObject

Object which contains information of the intermediate FITS file (to be written if intermediate = 1)

gsl_vector *record

GSL vector with signal values of input record

PulsesCollection *foundPulses

Input/output structure where the info about found pulses is stored

int pulseGrading(ReconstructInitSIRENA *reconstruct_init, int grade1, int grade2, int OFlength_strategy, int *pulseGrade, long *OFlength)

Located in file: tasksSIRENA.cpp

This function provides the pulse grade (Pileup=-2, Rejected=-1, HighRes=1, MidRes=2, LimRes=3, LowRes=4) and the optimal filter length by taking into account the info read from the XML file and the OFStrategy (FREE, BASE2, BYGRADE or FIXED).

Members/Variables

ReconstructInitSIRENA **reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values). In particular, this function uses the OFLength and grading

int grade1

Pulse duration (length of optimal filter applied)

int grade2

Difference between the start time of the pulse and the start time of the previous pulse

int OFlength_strategy

Same as OFStrategy (input)

int *pulseGrade

Pulse grade (output)

long *OFlength

Optimal filter length (= OFLength only if OFStrategy = FIXED and OFLength <= grade1) (output)

int readAddOrderParams(ReconstructInitSIRENA *reconstruct_init, fitsfile **inLibObject, double samprate, int eventcntLib, double estenergy, gsl_vector *pulsetemplate, gsl_matrix *covariance, gsl_matrix *weight, gsl_vector *pulsetemplateMaxLengthFixedFilter)

Located in file: tasksSIRENA.cpp

This function reads the library data, add new data and sort the data according to an energy-ascending order.

  • Declare variables
  • Load values already in the library
  • Add new values
  • Realign
  • Add intermeadiate values
  • Recalculate intermediate values of some new pairs
  • Write values in the library
  • Free allocated GSL vectors

Members/Variables

ReconstructInitSIRENA **reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values). In particular, this function uses FilterMethod, PulseLength, library_collection, monoenergy, mode and noise_spectrum

fitsfile **inLibObject

FITS object containing information of the library FITS file

double samprate

Sampling rate

int eventcntLib

Number of templates in the library

double estenergy

Pulse height of the template whose energy is going to be added to the library

gsl_vector *pulsetemplate

GSL vector with the pulse template whose energy is going to be added to the library

gsl_matrix *covariance

GSL matrix with covariance matrix of the energy which is going to be added to the library

gsl_matrix *weight

GSL matrix with weight matrix of the energy which is going to be added to the library

gsl_vector *pulsetemplateMaxLengthFixedFilter

GSL vector with the largeFilter-length template whose energy is going to be added to the library

int readFitsComplex(IOData obj, gsl_matrix **result)

Located in file: inoututils.cpp

This function reads values of a complex column of a FITS file. After that, the function puts them into a GSL matrix for an easier processing.

Members/Variables

IOData obj

Input object for complex FITS column

gsl_matrix **result

Output GSL matrix

int readFitsSimple(IOData obj, gsl_vector **result)

Located in file: inoututils.cpp

This function reads values of a simple column of a FITS file. After that, the function puts them into a GSL vector for an easier processing.

Members/Variables

IOData obj

Input object for simple FITS column

gsl_vector **result

Output GSL vector

extern_C_void reconstructRecordSIRENA(TesRecord *record, TesEventList *event_list, ReconstructInitSIRENA *reconstruct_init, int lastRecord, int nRecord, PulsesCollection **pulsesAll, OptimalFilterSIRENA **optimalFilter, int *const status)

Located in file: integraSIRENA.cpp

This function is the main wrapper function to detect, grade and calculate energy of pulses in input records.

  • Inititalize PulsesCollection structure

  • Check consistency of some input parameters

  • Detect pulses in input record (runDetect()). Save information of detected pulses

    • If PCA, pulses energies are already written in the pulsesAll structure
  • If in RECONSTRUCTION (mode = 1):

    • Filter record and calculate energy of pulses (runEnergy())
  • Populate output event list with pulses energies, arrival time and grading

Members/Variables

TesRecord *record

Instance of TesRecord structure that contains the input record

TesEventList *event_list

Instance of TesEventList structure that contains the information of the reconstructed pulses

ReconstructInitSIRENA *reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values)

int lastRecord

If record being analyzed is the last one, lastRecord = 1. Otherwise it is equal to 0

int nRecord

Input record number

PulsesCollection **pulsesAll

Member of PulsesCollection structure to successively store all the pulses used to create the library. Re-populated after each processed record.

OptimalFilterSIRENA **optimalFilter

Optimal filters used in reconstruction

int *const status

Input/output status

int RS_filter(gsl_vector *vector, double lrs, double lb, double B, double *pulseheight)

Located in file: pulseprocess.cpp

This function uses the running sum filter to find the pulse height. It always works in time domain.

A running sum filter, RS, is the sum of lrs digitized data samples. It is continuously updated upon the arrival of new data point. Simultaneously a baseline filter, B, is the sum of lb digitized data samples without pulses. The algorithm looks for the time when RS/lrs reaches its maximum. At that time RS is stored, RS_{max}, and the baseline is scaled with lrs, Bp ( Bp=B \cdot lrs/lb). Then, the pulse height related to the pulse pseudoenergy is given by:

Pulse height=\frac{RS_{max}-B_p}{lrs}

Members/Variables

gsl_vector *vector

Not filtered pulse (extracted from the record in getPulseHeight())

double lrs

Running sum length (samples)

double lb

Baseline averaging length (samples)

double B

In general, sum of the lb digitized data samples of a pulse-free interval immediately before the current pulse

double *pulseheight

Pulseheight of the input pulse

void runDetect(TesRecord *record, int lastRecord, PulsesCollection *pulsesAll, ReconstructInitSIRENA **reconstruct_init, PulsesCollection **pulsesInRecord)

Located in file: tasksSIRENA.cpp

This function is responsible for the detection in SIRENA, record by record. It is used both for library creation (mode = 0) and energy reconstruction (mode = 1) runnings.

Conditions:

Steps:

  1. Create library file if it is necessary: last record and mode = 0 (run createLibrary())

  2. Create intermediate output FITS file if required (createDetectFile())

  3. (Filter and) differentiate the models of the library (only for the first record in mode = 1 run (filderLibrary())

  4. Store the input record in invector (loadRecord())

  5. Convert I into R if EnergyMethod = I2R, I2RALL, I2RNOL or I2RFITTED (convertI2R())

  6. Process each record (proceRecord()):

    • (Low-pass filter and) differentiate
    • Find pulses
    • Load the found pulses data in the input/output foundPulses structure
    • Write test info in intermediate output FITS file if reconstruct_init->intermediate = 1 (writeTestInfo())
    • Write pulses info in intermediate output FITS file if reconstruct_init->intermediate = 1 (writePulses())

From this point forward, I2R, I2RALL, I2RNOL and I2RFITTED are completely equivalent to OPTFILT

  1. If last record in mode = 0 run:

  2. If last record and PCA:

    • In order to not have restrictions when providing (*reconstruct_init)->energyPCAx
    • Covariance data
    • Eigenvalues and eigenvectors
    • RSxN (S=2)
    • AE straight line: Pto0(x,y) and Pto10(x,y)
    • Calculus of the rotation angle
    • Rotation
    • Histograms of the two clusters (two energies)
    • Conversion factor from arbitrary unit to eV
    • Energy calculation
  3. Close intermediate output FITS file if it is necessary

Members/Variables

TesRecord *record

Member of TesRecord structure that contains the input record

int lastRecord

Integer to verify whether record is the last one (=1) to be read (and thus if library file will be created)

PulsesCollection *pulsesAll

Member of PulsesCollection structure to successively store all the pulses used to create the library. Re-populated after each processed record

ReconstructInitSIRENA **reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values)

PulsesCollection **pulsesInRecord

Member of PulsesCollection structure to store all the pulses found in the input record

void runEnergy(TesRecord *record, ReconstructInitSIRENA **reconstruct_init, PulsesCollection **pulsesInRecord, OptimalFilterSIRENA **optimalFilter)

Located in file: tasksSIRENA.cpp

This function calculates the pulse energy applying different methods (from EnergyMethod and OFNoise). It only runs in RECONSTRUCTION mode (mode = 1) (except to EnergyMethod = PCA).

  • Declare variables

  • Store the record in invector (loadRecord())

  • Subtract the baseline if EnergyMethod = OPTFILT and runF0orB0val = 1 (FilterMethod = B0)

  • Subtract the baseline if EnergyMethod = WEIGHT

  • Check Quality

  • For each pulse:

    • Establish the pulse grade (HighRes=1, MidRes=2, LimRes=3, LowRes=4, Rejected=-1, Pileup=-2) and the optimal filter length

    • Pulse: Load the proper piece of the record in pulse

    • If OFIter = 1, in the first iteration ( numiteration = 0) the values of maxDER and maxDERs are used in find_matchedfilter(), find_matchedfilterDAB(), find_optimalfilter(), find_optimalfilterDAB() or find_Esboundary() getting the values of the energies which straddle the maxDER (Ealpha and Ebeta). It will have more iterations if the calculated energy is out of [Ealpha, Ebeta]. If energy is in [Ealpha, Ebeta] the iterative process stops.

      • If EnergyMethod = OPTFILT (or I2R, I2RALL, I2RNOL, I2RFITTED) and OFLib = 0 and OFNoise = NSD:

      • If EnergyMethod = OPTFILT (or I2R, I2RALL, I2RNOL, I2RFITTED) and OFLib = 1 and OFNoise = NSD:

      • If EnergyMethod = WEIGHT or WEIGHTN:

        • Get the indexes of the two energies which straddle the pulse

        • If EnergyMethod = WEIGHTN and OFLib = 1:

          • Choose the base-2 system value closest (lower than or equal) to the pulse length
          • find_prclwn() to find the appropriate values of the PRECALWN HDU (PRCLx columns)
      • If EnergyMethod = OPTFILT (or I2R, I2RALL, I2RNOL, I2RFITTED) and OFLib = 1 and OFNoise = WEIGHTM:

        • Choose the base-2 system value closest (lower than or equal) to the pulse length
        • find_prclofwm() to find the appropriate values of the PRCLOFWM HDU (OFWx columns)
      • Calculate the energy of each pulse

      • If using lags, it is necessary to modify the tstart of the pulse and the length of the filter used

    • In order to subtract the pulse model, it has to be located in the tstart with jitter and know its values in the digitized samples

    • Subtract the pulse model from the record

    • Write info of the pulse in the output intemediate file if intermediate = 1

  • Free allocated GSL vectors

Members/Variables

TesRecord **record

Structure that contains the input ADC record

ReconstructInitSIRENA **reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values)

PulsesCollection *pulsesInRecord

Collection of pulses found in the current record

OptimalFilterSIRENA **optimalFilter

Optimal filters used in reconstruction

int shiftm(gsl_vector *vectorin, gsl_vector *vectorout, int m)

Located in file: tasksSIRENA.cpp

This function returns as vectorout the vectorin delayed m samples.

Members/Variables

int m

Delay in samples

gsl_vector *vectorin

GSL vector with input vector

gsl_vector *vectorout

GSL with input vector (vectorin) delayed m samples

int shift_m(gsl_vector *vectorin, gsl_vector *vectorout, int m)

Located in file: tasksSIRENA.cpp

This function returns as vectorout the vectorin moved forward by m samples.

Members/Variables

int m

Advance in samples

gsl_vector *vectorin

GSL vector with input vector

gsl_vector *vectorout

GSL with input vector (vectorin) moved forward m samples

int toGslMatrix(void **buffer, gsl_matrix **matrix, long numCol, int numRow, int type, int eventini)

Located in file: inoututils.cpp

The function puts the values of the input buffer into an output GSL matrix. Columns and rows are input parameters.

Members/Variables

void **buffer

Input buffer with data

gsl_matrix **matrix

Output GSL matrix

long numCol

Number of columns

int numRow

Number of rows

int type

FITS type (TINT, TSHORT, TDOUBLE, etc.)

int eventini

Initial event to start writing

int toGslVector(void **buffer, gsl_vector **array, long nevent, int eventini, int type)

Located in file: inoututils.cpp

The function puts the values of the input buffer into an output GSL vector.

Members/Variables

void **buffer

Input buffer with data

gsl_vector **array

Output GSL vector

long nevent

Number of elements to store

int eventini

Initial element number

int type

FITS type (TINT, TSHORT, TDOUBLE, etc.)

int vector2matrix(gsl_vector *vectorin, gsl_matrix **matrixout)

Located in file: tasksSIRENA.cpp

This function converts an input n^2 vector into an output [n \times n] square matrix. It puts the first n elements of the vector in the first row of the matrix, the second group of n elements (from n to 2n-1) of the vector in the second row and so on.

Members/Variables

gsl_vector **vectorin

GSL input vector whose length is n^2

gsl_matrix *matrixout

GSL output square matrix [n \times n]

int weightMatrix(ReconstructInitSIRENA *reconstruct_init, bool saturatedPulses, PulsesCollection *pulsesAll, PulsesCollection *pulsesInRecord, long nonpileupPulses, gsl_vector *nonpileup, gsl_vector *pulseaverage, gsl_matrix **covariance, gsl_matrix **weight)

Located in file: tasksSIRENA.cpp

[Fixsen+04]

This function calculates the weight matrix by using the non piled-up pulses found in all the records, stored in pulsesAll (previous records) and pulsesInRecord (current record). The weight matrix of each energy (and other intermediate values) will be stored in the library by the function fillInLibraryData().

Definitions:

S_i^p: Value of the ith-sample of the pulse number p

M_i^p: Value of the ith-sample of the model number p (model= pulseaverage):

M_i = <S_i> = (1/N)\sum_{p=1}^{N}S_i^p

N: number of non piled-up pulses

& D_i = S_i - M_i \\
& V_{ij} = <D_iD_j> = E[(S_i-M_i)(S_j-M_j)] = (1/N)\sum_{p=1}^{N}(S_i^p-M_i^p)(S_j^p-M_j^p) \\
& V = \left[\begin{matrix} <D_1D_1> & <D_1D_2> & ... & <D_1D_n> \\
<D_2D_1> & <D_2D_2> & ... & <D_2D_n> \\
....  &  ....  & ... &  ....  \\
<D_nD_1> & <D_nD_2> & ... & <D_nD_n>\end{matrix}\right]

where n is the PulseLength and thus V = [n \times n].

The weight matrix W = [V]^{-1}.

Steps:

  • Calculate the elements of the diagonal of the covariance matrix

  • Calculate the elements out of the diagonal of the covariance matrix

  • If saturated pulses \Rightarrow Covariance matrix is a singular matrix \Rightarrow Non invertible

    In order to allow the covariance matrix to be inverted \Rightarrow Replacing 0’s (0’s are due to the saturated values, equal in the pulse and in the model)

    • Elements of the diagonal: Generating a random double f_1 between a range (fMin,fMax) (-NoiseStd,NoiseStd) to replace 0’s with f_1^2
    • Elements out of the diagonal: Generating two random doubles f_1 and f_2 between a range (fMin,fMax) (-NoiseStd,NoiseStd) to replace 0’s with f_1 \cdot f_2
  • Calculate the weight matrix

Members/Variables

ReconstructInitSIRENA **reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values). In particular, this function uses PulseLength and EnergyMethod

bool saturatedPulses

If true, all the pulses ( mode = 0 \Rightarrow all the pulses have the same energy) are saturated

PulsesCollection *pulsesAll

Collection of pulses found in the previous records

PulsesCollection *pulsesInRecord

Collection of pulses found in the current record

long nonpileupPulses

Number of non piled-up pulses

gsl_vector **pulseaverage

GSL vector with the pulseaverage (= template = model) of the non piled-up pulses

gsl_matrix **covariance

GSL matrix with covariance matrix

gsl_matrix **weight

GSL matrix with weight matrix

int writeFilterHDU(ReconstructInitSIRENA **reconstruct_init, int pulse_index, double energy, gsl_vector *optimalfilter, gsl_vector *optimalfilter_f, gsl_vector *optimalfilter_FFT, fitsfile **dtcObject)

Located in file: tasksSIRENA.cpp

This function runs in RECONSTRUCTION mode and writes the optimal filter info (in the FILTER HDU) for each pulse if intermediate = 1.

  • Declare variables
  • Open intermediate FITS file
  • Create the FILTER HDU if it is the first pulse
  • Write data
    • OPTIMALF column
    • OFLENGTH column
    • FREQ column
    • OPTIMALFF column
    • ENERGY column
  • Modify the Primary HDU
  • Close intermediate output FITS file if it is necessary
  • Free memory

Members/Variables

ReconstructInitSIRENA **reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values). In particular, this function uses detectFile, clobber and PulseLength

int pulse_index

Index of the pulse whose info is going to be written (to know if it is the first pulse)

double energy

Estimated energy (eV)

gsl_vector *optimalfilter

Optimal filter in time domain

gsl_vector *optimalfilter_f

Frequency axis of the optimal filter spectrum

gsl_vector *optimalfilter_FFT

Optimal filter spectrum

fitsfile **dtcObject

Fitsfile object for intermeadiate file name

int writeFitsComplex(IOData obj, gsl_matrix *matrix)

Located in file: inoututils.cpp

This function reads values of a GSL matrix. After that, the function puts them into a complex column of the output FITS file.

Members/Variables

IOData obj

Object for FITS column to be written

gsl_matrix *matrix

Input GSL matrix with data

int writeFitsSimple(IOData obj, gsl_vector *vector)

Located in file: inoututils.cpp

This function reads values of a GSL vector. After that, the function puts them into a column of the output FITS file.

Members/Variables

IOData obj

Object for FITS column to be written

gsl_vector *vector

Input GSL vector with data

int writeLibrary(ReconstructInitSIRENA *reconstruct_init, double samprate, double estenergy, gsl_vector *pulsetemplate, gsl_matrix *covariance, gsl_matrix *weight, bool appendToLibrary, fitsfile **inLibObject, gsl_vector *pulsetemplateMaxLengthFixedFilter)

Located in file: tasksSIRENA.cpp

This function writes the library (reordering if it is necesary and calculating some intermediate parameters)

  • Adding a new row to the library
  • Write the first row of the library

Members/Variables

ReconstructInitSIRENA **reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values)

double samprate

Sampling rate

double estenergy

Pulse height of the template whose energy is going to be added to the library

gsl_vector *pulsetemplate

GSL vector with the pulse template whose energy is going to be added to the library

gsl_matrix **covariance

GSL matrix with covariance matrix

gsl_matrix **weight

GSL matrix with weight matrix

bool appendToLibrary

true if adding a new row to the library and false if it is the first row to be added

fitsfile **inLibObject

FITS object containing information of the library FITS file

gsl_vector *pulsetemplateMaxLengthFixedFilter

GSL vector with the largeFilter-length pulse template whose energy is going to be added to the library

void writeLog(FILE *fileRef, string type, int verbosity, string message)

Located in file: inoututils.cpp

This function includes the processing of the each level of message in the log file and the output screen:

  • Verbosity = 0 \Rightarrow The log file and the output screen include Errors
  • Verbosity = 1 \Rightarrow The log file and the output screen include Errors and Warnings
  • Verbosity = 2 \Rightarrow The log file and the output screen include Errors, Warnings and Alerts
  • Verbosity = 3 \Rightarrow The log file and the output screen include Errors, Warnings, Alerts and Log messages

Members/Variables

FILE *fileRef

File reference to log file

string type

String to indicate error type “Error”, “Warning”, “Alert”,”Log” or “OK”

int verbosity

Integer value for verbosity

string message

String message to print

int writePulses(ReconstructInitSIRENA **reconstruct_init, double samprate, double initialtime, gsl_vector *invectorNOTFIL, int numPulsesRecord, gsl_vector *tstart, gsl_vector *tend, gsl_vector *quality, gsl_vector *taurise, gsl_vector *taufall, fitsfile *dtcObject)

Located in file: tasksSIRENA.cpp

This function writes the data of the pulses found in the record in the intermediate FITS file (in the PULSES HDU). The pulses info given is: TSTART, I0 (the pulse itself), TEND, TAURISE, TAUFALL and QUALITY.

Members/Variables

ReconstructInitSIRENA **reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values). In particular, this function uses detectFile, PulseLength and clobber.

double samprate

Sampling rate (to convert samples to seconds)

double initialtime

Starting time of the record (in order to calculate absolute times)

gsl_vector *invectorNOTFIL

GSL vector with the original record (neither low-pass filtered nor differentiated)

int numPulsesRecord

Number of pulses found in the record

gsl_vector *tstart

GSL vector with the start times of the found pulses

gsl_vector *tend

GSL vector with the end times of the found pulses

gsl_vector *quality

GSL vector with the quality of the found pulses

0 \Rightarrow Standard (good) pulses

1 \Rightarrow Truncated pulses

10 \Rightarrow Saturated pulses

11 \Rightarrow Truncated and saturated pulses

gsl_vector *taurise

GSL vector with the rise time constants of the found pulses (to be done)

gsl_vector *taufall

GSL vector with the fall time constants of the found pulses (to be done)

fitsfile *dtcObject

Object which contains information of the intermediate FITS file

int writeTestInfo(ReconstructInitSIRENA *reconstruct_init, gsl_vector *recordDERIVATIVE, double threshold, fitsfile *dtcObject)

Located in file: tasksSIRENA.cpp

This function writes the TESTINFO HDU in the intermediate FITS file. The written columns are FILDER (low-pass filtered and differentiated record) and THRESHOLD.

Members/Variables

ReconstructInitSIRENA **reconstruct_init

Member of ReconstructInitSIRENA structure to initialize the reconstruction parameters (pointer and values). In particular, this function uses detectFile.

gsl_vector *recordDERIVATIVE

GSL vector with input record (low-pass filtered and) differentiated

double threshold

Threshold value used to find pulses

fitsfile dtcObject

Object which contains information of the intermediate FITS file