Home

BuL Virgo Gravitational Wave Searches Library for Burst

1. amplitude frequency Figure 2 8 Amplitude and phase of the transfer function of a low pass band filter using the Butterworth approximation N 20 bilir itr 2000 4000 8000 10000 12000 14000 16000 18000 20000 10 10 frequency Figure 2 9 Data have been generated following the Virgo sensitivity curve black curves and are then low pass filtered The red curves represent the filtered data upper plot whose PSD is represented below 48 2 5 4 2 Butterworth high pass filter order 30 tamalitude 10 E q i 34 E eth a E ee un a 19 E ojear cid LE da a E BAA A 10 10 frequency 3 1 1 2 CE re i 10 10 of Figure 2 10 Amplitude and phase of the transfer function of a high pass band filter using the Butterworth approximation N 30 1 wild il l A Ah W Y poetiiitiirirtiiitrrrtipitiii 2000 4000 6000 8000 10000 12000 14000 16000 18000 2000 6000 18000 20000 2 10 10 frequency Figure 2 11 Data have been generated following the Virgo sensitivity curve black curves and are then high pass filtered The red curves represent the filtered data upper plot whose PSD is represented below 49 2 5 4 3 Butterworth pass band filter 10 10 frequency
2. y 10 r mo a Ey ey eee ooo RS 10 9 amplitude ell dal pepe pa SR A Na peepee ry ES A a TITTTTTTTTTTTTTTTTTT TT TT o o o 10 10 10 t o Figure 2 6 Transfer function amplitude and phase of a multi pole and zero filter designed using the IIR class The x axis represents frequency in Hz 41 2 5 3 Class IIR Include file LinearFilter IIR hxx Constructors IIR Vector lt double gt amp aZerosFrequency Vector lt double gt amp aZerosQ Vector lt double gt amp aPolesFrequency Vector lt double gt amp aPolesQ double aGain double aSampling Destructors IIR O Public methods void SetGain double aValue int SetGain const double aValue const int aRampTime double GetCurrentGain void SetGainAtFrequency double aValue double aFrequency double GetGainAtFrequency double aFrequency void ResetHistory void SetHistory int aDepth Vector lt double gt amp aX Vector lt double gt amp aY void GetHistory int aDepth Vector lt double gt amp aX Vector lt double gt amp aY void PrintCoefficients double Filter double alnput int Filter Vector lt double gt amp alnput Vector lt double gt amp a0utput int GetTransferFunction Vector lt double gt amp aModulus Vector lt double gt amp aPhase double aDeltaFrequency int GetTransferFunction Vector lt complex lt double gt gt amp aTransfer
3. and GetSamplingFrequency methods allow to change and get the number of poles of the AR filter 85 Program example include Noise_ARMA hxx include include Spectral_Spectral hxx include ITF_Virgo hxx using namespace _BuL int main int arg char argc int i j sampling size int poles_nb zeros_nb double f_cut sampling 20000 poles_nb 38 zeros_nb 16 size sampling 2 f_cut 10 double factor 1 e40 Virgo PSD Vector lt double gt PSD 100 size Virgo itf sampling itf GetOneSidedPSD PSD DataContainer_Vector hxx precision de 0 01 Hz ARMA model sampling poles_nb zeros_nb model ComputeARMA PSD factor f_cut model SaveARMA test dat model ComputeARMA ITF Virgo f_cut model SaveARMA test2 dat return 0 86 2 7 6 Class GaussianLine 2 7 6 1 Class Line Include file Noise Line hxx Constructors Line Line bool aDebug Line int aSampling double aFrequency double aQ double aMass double aTemperature double aCoefficient Line int aSampling double aFrequency double aSigmaProcess double aCoefficient Line Line amp aLine Destructors Line O Public methods int SetThermal int aSamplingFrequency double aFrequency double al double aMass double aTemperature double aCoefficient int SetThermal double aFrequency double aQ double aMass double aTemperature double aCoefficient int SetElectr
4. double GetSNR double GetHrss double GetDistance int SetSNR double aSNR int SetHrss double aHrss int SetDistance int aDistance double GetSNRForHrss double aHrss double GetHrssForSNR double aSNR int SetPSD int aSampling Vector lt double gt amp aPSD int SetITF ITF amp alTF int SetSigmaNoise int aSampling double aSigma double GetSigmaNoise int GetSamplingFrequency int GetData double aPosition Vector lt double gt amp aData void SetDebug bool aValue bool GetDebug 101 2 9 BuT burst template generators 2 9 1 Class TGaussian This class provides methods to define and generate templates of Gaussian form in time space as well as in Fourier space The Gaussian template wave form is given by g t e 2 27 where to is chosen as half of the window size in such a way the template is centered By default the Gaussian template is not normalized since it depends on a given PSD Include file BuT_TGaussian hxx Constructors TGaussian double aSigmaMin double aSigmaMax double aEpsilon int aSampling TGaussian double aSigmaMin double aSigmaMax double aEpsilon int aSampling bool aPadding Public methods int GetNTemplate double GetSigmaTemplate int alndex Vector lt double gt amp GetSigmaTemplate int alndex void SetDebug bool aValue bool GetDebug int GetWindowSize void SetWindowSize int aValue int GenerateWaveForm int GenerateFourier WaveForm Matrix l
5. complex lt double gt ali 0 bb i complex lt double gt b i 0 alsizeti 0 zero padding b sizeti 0 zero padding aa size i complex lt double gt 0 0 zero padding bb size i complex lt double gt 0 0 zero padding cout lt lt Correl WITH ZERO PADDING lt lt endl Correl toto size true toto Correlation a b c toto Correlation aa bb cc toto AutoCorrelation a c toto AutoCorrelation aa cc toto Convolution a b c toto Convolution aa bb cc return 0 14 2 2 3 Class GPSTime This class provides some methods to convert a GPS time into UTC time or local system time and vice versa When converting a GPS Time into a UTC time one has to take into account the leap seconds which are regularly added to UTC time in order to ensure that the difference between a uniform time scale defined by atomic clocks TAT time does not differ from the Earth s rotational time by more than 0 9 second 4 The UTC epoch is January 1 1970 The GPS epoch is January 6 1980 3600 24 365 10 2 5 The first leap second added to UTC was introduced on June 30 1972 The last leap second added to UTC was introlduced on December 31 1999 As of January 1 1999 TAI is ahead of UTC by 32 seconds TAI is ahead of GPS by 19 seconds GPS is ahead of UTC by 13 seconds All these information are stored in a file S BUMROOT data leapseconds txt which is read when calling the GP
6. 127 2 11 4 Class Example or how to code a new Burst algorithm This class is a skeleton class derived from Detector which can be taken as the starting point to define anew Burst algorithm class All the virtual methods of the class Detector are implemented but are empty The test application source file is test _BuF_Example cxx Include file Constructors Public methods BuF_Example hxx Example string aName bool Configure string aFileName bool Amplitude Vector lt double gt amp aData Matrix lt double gt amp aAmplitude int Threshold Matrix lt double gt amp aAmplitude list lt Event gt amp aEventList bool Detect Vector lt double gt amp aData Matrix lt double gt amp aAmplitude list lt Event gt amp aEvent List int End int GetNWindow int GetWindow int alndex void SetDebug bool avalue bool GetDebug The method GetNWindow returns the number of window which will be used in paralell The method GetWindow int alndex returns the size of the window which corresponds to the nimber alndez test_BuF_Example cxx program include lt math h gt include BuF_Example hxx tinclude DataContainer_Matrix hxx tinclude Noise_Gaussian_MT hxx include BuS_Gaussian hxx using namespace _BuL int main int narg char argc int i int size 2048 Matrix lt double gt Amplitude 10 size Vector lt double gt h size Vector lt double gt g si
7. 2 81448e 07 s 2 1 0 00159155 s 0 000253303 s 2 gt Filter No 1 Filter in Laplace space 1 0s 0 s72 1 0 000159155 s 2 53303e 06 s 2 gt Filter No 2 Filter in Laplace space 1 O0s 0s72 1 1 59155e 07 s 2 53303e 08 s 2 hhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhh hh The Filter methods perform the filtering transform either on a given input sequence value or on a full sequence vector Note that in the case of a full sequence vector if one changes the filter gain it will be effective at the begining of the next sequence The GetTransferFunction method computes the amplitude and the phase of the transfer function using Eq 2 5 One has to specify the frequency step to fill out the vectors there exists also the method which provides the complex transfer function The LowPassToHighPass Transform method transforms a low high pass filter into a low high pass filter as explained in 2 5 4 Example Low pass filter using IIR class include lt stdio h gt include lt stdlib h gt include lt math h gt include lt fstream gt include DCPP_Vector hxx include LinearFilter_IIR hxx include Spectral_Taper hxx include Spectral_Spectral hxx using namespace _BuL low pass filter 43 int main int argc char argv ofstream file IIR filter double output int i double sampling 20000 double gain int datasize 20000 Vector lt double gt module Vect
8. T NF D y2n 1 Note that the output for each window is normalized and the threshold level may be the same for all the windows The test application source file is test_BuF_NF cxx Include file BuF NF hxx Constructors NF string aName Public methods int Configure string aFileName int Configure Vector lt int gt amp aWindowSize Vector lt double gt amp aThreshold int aSampling int Amplitude Vector lt double gt amp aData Matrix lt double gt amp aAmplitude int Threshold Matrix lt double gt amp aAmplitude list lt Event gt amp aEventList int Threshold Matrix lt double gt amp aAmplitude double amp aGPsS list lt Event gt amp aEventList int Detect Vector lt double gt amp aData Matrix lt double gt amp aAmplitude list lt Event gt amp aEventList int Detect Vector lt double gt amp aData double amp aGPS list lt Event gt amp aEventList int End int SetThreshold Vector lt float gt aThreshold int GetThreshold Vector lt float gt amp aThreshold int GetNWindow int GetWindow int alndex void SetDebug bool avalue bool GetDebug int NormFilter Vector lt double gt amp aData int aWindow Vector lt double gt amp aOutput The method SetThreshold changes the initial values of the thresholds defined for each window The method GetThreshold returns the value of the thresholds The method GetNWindow returns the number of windows which are used in
9. class Include file TTF_Virgo hxx int aSamplingFrequency ITF Detector aName ITF Detector aName int aSamplingFrequency Constructors Virgo Virgo Virgo Virgo PS Can i aan Public methods int GetTwoSidedPSD Vector lt double gt amp aData int GetOneSidedPSD Vector lt double gt amp aData int GetTwoSidedASD Vector lt double gt z aData int GetOneSidedASD Vector lt double gt amp aData int GetOneSidedASD Vector lt double gt amp aData double aFactor int GetT woSidedASD Vector lt double gt amp aData double aFactor int GetOneSidedPSD Vector lt double gt amp aData double aFactor int GetT woSidedPSD Vector lt double gt amp aData double aFactor void SetSamplingFrequency int aValue int GetSamplingFrequency double GetLatitude double GetLongitude double Get Azimuthal double GetArmAngle void SetLatitude double aValue void SetLongitude double aValue void SetAzimuthal double aValue void SetArmAngle double aValue void SetDebug bool aValue bool GetDebug The Virgo constructor consider the ITF Virgo detector with its nominal fo 20 kHz sampling frequency The three other constructors allow to choose another Virgo detector and or set up a sampling frequency if needed otherwise the default value is 20 kHz The methods are described in the section of the ITF class Program exemple include ITF_Virgo hxx using namespace _BuL int main int arg char argc int i j
10. just to check the normalization filter SetTemplate patron just to check the PSD setting filter SetPSD psd_twosided autocorrelation data patron filter Filter data output output WriteFile match dat return 0 58 2 5 8 Class Goertzel The Goertzel algorithm computes the DFT coefficient associated to the desired frequency without the frequency sampling restriction inherent to the DFT w 2akfo N The computation is propertionel to N let s recall that FFT algorithm computation goes as NlogN The main advantage of the Goertzel algorithm is to restrict the computation to the desired frequencies Given an input data size N that still constraints the frequency resolution the sampling frequency fo and the frequency of interest f the algorithm computes the DFT coefficient at index k N f fo The Goertzel algorithm can be viewed as the response of a system to a finite length input x n sequence 9 That s the reason why the algorithm is commonly implemented with a second order ITR bandpass filter with the following transfer function 1 71 H z 2 20 2 1 2cos ZE z 1 272 29 In the time domain the DFT kt coefficient X k is given by X k s N e2 Ns N 1 the sought DFT coefficient 2 21 27k spin z n 2 cos skin 1 s n 2 recursion relation 2 22 sp 2 0 initial conditions 2 23 sk 1 0 2 24 Include file LinearFilter_Goertzel hxx Const
11. Vector lt complex lt double gt gt amp aDataT wo Vector lt complex lt double gt gt amp aData int AutoCorrelation Vector lt double gt amp aDataOne Vector lt double gt amp aData int AutoCorrelation Vector lt complex lt double gt gt amp aDataOne Vector lt complex lt double gt gt amp aData int Convolution Vector lt double gt amp aDataOne Vector lt double gt zaData Two Vector lt double gt amp aData int Convolution Vector lt complex lt double gt gt amp aDataOne Vector lt complex lt double gt gt amp aDataTwo Vector lt complex lt double gt gt amp aData Integral gt ay aDatali x aStep N aData size IntegralNorm Ea aData i aStep N aData size IntegralModule Ni laDatali x aStep N aData size GetThresholdKhi2 returns the value of the threshold corresponding to the given false alarm rate aFalseAlarm asuming a x distribution of aN degrees of freedom GetThresholdGaussian returns the value of the threshold corresponding to the given false alarm rate aFalseAlarm asuming a gaussian distribution of width aSigma The folowing functions have been extracted from Numerical Recipes double Erfc double aX erfc x double Gammln double aX InIT x double Gammaq double aA double aX Q a x a gt 0 double Gammp double aA double aX P a x a gt 0 void Gcf double amp aGammef double aA double aX double amp aGln Q a x a gt 0 x gt a 1 void Gser double amp aGamser
12. const string aFileName int Add const double aValue const double aWeight The class is derived from the Histogram class to provide an alternative data accumulation Instead of counting the number of events falling in a range of values events are segregated following a first variable while a second one is saved so as to obtain its mean value and variance The Define method calls the Hixstogram class method but does not require a weight value as it is this weight that is modified for each event Add aggregates events according to the first variable with the second variable stored to get mean and WriteSigmaHistogramFile returns the variance for each histogram bin using the Histogram file format The mean values are obtained using WriteHistogramFile inherited from Histogram 21 2 2 7 DataCard The DataCard class may be used to read parameters stored in a file Each line of the file starting with the keyword will be considered These lines contain a set of values associated to a keyword as shown in the following data card example All other starting lines are seen as comments lines Include file BuM_DataCard hxx Constructors DataCard string aFileName Destructors DataCard Public methods int ReadCard int GetParam string amp aKey int amp alndex int amp aValue int GetParam string amp aKey int amp alndex float amp aValue int GetParam string aKey int alndex double amp aValue int GetParam s
13. dc GetParam ALGO O name cout lt lt name lt lt endl dc GetParam datasize 0 test dc GetParam TYTO 0 test cout lt lt write out dataCard lt lt endl dc WriteCard return 0 3 23 2 3 Spectral This package contains several classes dedicated to non parametric spectral estimation algorithms Most of the content of this package are well explained in 5 The power spectral density function of a process x t is defined as the Fourier Transform of the autocor relation function A T of z t S f A PSD f Fa AQ f ERR X The amplitude spectral density function is defined as the square root of the power spectrum 2 3 1 Non parametric spectral estimation One usually tries to estimate the Power Spectral Density Function S f of a physical process using a data set composed of a given number N of discrete value xy t This induces a bias in all the PSD estimators which can be attributed to the existence of sidelobs in the Fej r s kernel F which appears in the simple PSD estimator periodogram EN EE a eS gt F t ESA where _ sin Naf Nsin2af The true PSD is convolved with the Fej r s kernel Fourier transform of the square function which is introduced by the size limitation of the data set Note that the sidelobs amplitude of the Fej r s kernel compared to the principal lob decreases when N is increasing The sidelobs of the Fej r s kernel introduce a leak
14. f Le e Nis odd n x is the greater integer which is smaller than x index 0 FT for f 0 index 1 n DFT output for f 4 nf index n 1 N 1 DFT output for f n amp f k In the case of the Forward DFT of a real data vector the output can be either a complex data vector whose format has been described above or an half complex data vector according to the following scheme 4One has to note that in the case of a real data vector DFT the use of the half complex data format avoids a copy of the output vector although it is unavoidable when one requests a complex data vector e N is even index 0 FT for f 0 index 1 real part of the FT for f Lo Hii Ni index 1 N 1 imaginary part of the FT for Y 17 Il i e Nis odd n x is the greater integer which is smaller than x index 0 FT for f 0 index 1 n real part of the FT for f LL nf index n 1 N 1 imaginary part of the FT for nfe m Le Include file BuM_FFT hxx Constructors FFT int aSize FFT Type aType FFT int aSize FFT Type aType FFT Direction aDirection FFT int aSize FFT Type aType FFT Wisdom aWisdom FFT int aSize FFT Type aType FFT Direction aDirection FFT Wisdom aWisdom FFT int aSize FFT Type aType FFT Wisdom aWisdom char aFileName FFT int aSize FFT Type aType FFT Direction aDirection FFT Wisdom aWisdom char aFileName Pu
15. frequency Figure 2 12 Amplitude and phase of the transfer function of a pass pass band filter using the Butterworth approximation N 10 for the low pass filter and N 22 for the high pass filter x 10 di n alya frequency Figure 2 13 Data have been generated following the Virgo sensitivity curve black curves and are then band pass filtered The red curves represent the filtered data upper plot whose PSD is represented below 50 2 5 5 IIR filter using the Chebyshev approximation ol 2 5 6 Class Butteworth This class implements 3 kind of filters e low pass filter e high pass filter e pass band filter UN TF poss TF poss a TFs frequency Figure 2 14 Transfer function of a low pass filter left and a high pass filter right TF i gt T top low T pass Pe T stop high frequency Figure 2 15 Transfer function of a pass band filter A low pass and a high pass filter are defined by their stop band and pass band frequencies and the transfer function values corresponding to these frequencies see Fig 2 14 A pass band filter is defined by four frequencies and 3 transfer function values as given in Fig 2 15 52 Include file LinearFilter_Butterworth hxx Constructors Butterworth double aSampling double aFrequencyPass double aFrequencyStop double aTFPass double aTFStop Butterworth double aSampling double aFrequencyStopL
16. list lt Event gt macroEvent clusterBase eventList cout lt lt Macro events found lt lt endl cout lt lt macroEvent lt lt endl 125 produces Ts Te 2446 2469 Tmax 2456 Amax 7 18573 Template 1 T 2 GPS 0 Tolerance 40 Micro event template 3 5 Ts Te 2566 2625 Tmax 2571 Amax 6 25865 Template 1 T 3 GPS 0 Tolerance 40 Micro event template 0 1 4 Ts Te 2865 2865 Tmax 2865 Amax 5 05144 Template 0 T 1 GPS O Tolerance 10 Micro event template none Ts Te 11041 11041 Tmax 11041 Amax 5 48995 Template 0 T 1 GPS 0 Tolerance 10 Micro event template none The SortStart methods sort in increasing start time order a list of events The ClusterBase clusterizes a list of events produced by a basic filter according to the following procedure e sort in the increasing start time order the list of event The combine any ordered couple i j of events if Tipa of gt Thart the 2 events are then combined to define a new event k Teor T tart 5 end max T 4 Tina a Aras max Aj ozs AP ie TE corresponds to AF az oh max 0 01 Template Numero 1 Template Number 1 Micro event list contains the 2 events remove the events i and j from the list of events add the event k in the list of events continue to cluster the modified list of events ClusterCorrelatorTight clusterizes a list of events produced by a correlator filte
17. method computes for all the windows or templates the amplitude of the filter on the vector of data It returs an array windows templates number x data vector size containing the amplitude values The Threshold method compares the threshold level defined for each window or only once to the amplitude vector When the amplitude overcomes the threshold level that defines the starting time counted in bin of a micro event The ending time of this micro event is defined as the next count for which the amplitude decreases below the threshold value One can note that the previous definition of a starting time and an ending time for the micro event has very little physical meaning in the case of the correlator filters The physical width of the signal will rather be given by the width of the template which gives the highest output value in the case of several templates overcoming the threshold All the other fields see section 2 11 3 of the micro event are defined according to obvious definition except the tolerance width which depends on the filter familly e case of the basic filters NF MF ALF the tolerance width is defined as the size of the window which corresponds to this micro event e case of the Peak Correlator filter PC the tolerance width is defined as twice the sigma value of the corresponding gaussian template e case of the Damped Sine Correlator filter DSC the tolerance width is defined by the value o
18. the physical parameter space The difference between the ellipsoid and the real border of the efficiency area is not really taken into account in practical cases for one step searches MM is close to 1 and thus the ellipsoid is assumed to be a good approximation of the ambiguity surface Moreover the overlapping between close ellipses is an advantage to discard such a problem Clearly the value of MM from which the ellipsoid approximation of the efficiency area is no longer valid depends on the precise tiling problem considered 2 10 1 3 Tiler Algorithm It is well known that the optimal tiling of an infinite plane by identical disks of radius R consists in placing their centers on an hexagonal lattice separated by a distance of V3R see Figure 2 18 In this way the overlapping is minimized This property of circular tiling is used by the algorithm on the following way Once the location of a template ko in P has been chosen one computes its efficiency area Eg Through a simple plane transformation Ey becomes a circle Co of unit radius and the six neighbor centers are placed on the regular hexagonal lattice previously described This choice is nearly optimal if the shape of the efficiency areas is slowly varying with respect to their characteristic dimensions Using the inverse transformation the six new center positions are given in P Then each ellipse associated with a new center is tested in order to check if it must be kept The condit
19. which means that the loss of the SNR due to the mismatch between the template and the signal must be kept below 1 MM One can note that with this definition a of course pessimistic estimator of the fraction of false dismissals is 1 MM 2 31 For instance with MM 97 one has 10 This value of MM is usually found in the literature as the correspondence with this particular value of is easy This quantity M M is the only input parameter of the tiling procedure it allows one to define precisely the efficiency area E X MM of the template ky by the following equation 1 3 uv dM d lt 1 MM 2 32 The area of P including all the vectors ky i as which match this inequality is the inner part of an ellipsoid centered on X whose proper volume scales as 1 MM yer 2 The average fraction of recovered SNR for physical signals with parameters belonging to E X MM is at least equal to MM 107 Figure 2 18 Optimal tiling of an infinite plane by identical circles the centers belong to an hexagonal lattice 2 10 1 2 The tiling problem Tiling the parameter space consists of finding a set of N filters kx n so that the union of the 1 lt p lt surfaces EQ M M _ completely covers P The aim is to achieve this task by using as few 1 lt p lt templates as possible to keep the computing cost manageable The main problems encountered by any tiling procedure are overlapping gaps areas lost beyond
20. 0 The SetSeed allows to increment the internal seed number of 1 unit This method can be called at any instant The SetSeedLocalTime allows to increment the seed number used by drand48 of a given number which depends on the local time with a precision of 1 s The internal seed number is also increment of one unit This method can be called at any instant The GetRandomData returns a random variable whose value lays in the segment 0 1 Note that there is no change of a new sequence seed at each call of the function The GetRandomData aXmin aXmaz returns a random variable whose value lays in aX min aX maz Note that there is no change of a new sequence seed at each call of the function The GetRandomData Vector lt T gt aData fills a vector with random variable whose value lays in the segment 0 1 Note that contrary to the previous functions there is a change of a new sequence seed at each call of the function For that purpose the method SeedLocalTime is used 73 The GetGaussianData methods returns a Gaussian random variable whose value lays in the segment 0 1 The Gaussian variable is built from the random variable using the Box Muller algorithm described in 11 Note that there is no change of a new sequence seed at each call of the function 2 7 1 2 Class Random_MT Include file Noise Random MT hxx Constructors Random_MT Public methods void GetRandomData Vector lt double gt amp aData void
21. 095 M Beccaria E Cuoco G Curci Adaptive Identification of VIRGO like noise spectrum VIR NOT PIS 1390 096 and VIR NOT PIS 1390 095 T Damour B R Iyer B Sathyaprakash Phys Rev D63 044023 2001 T Zwerger E Mueller Astron Astrophys 320 209 1997 http www mpa garching mpg de Hydro GRAV index html H Dimmelmeier J A Font E Mueller Astron Astrophys 388 917 2002 H Dimmelmeier J A Font E Mueller Astron Astrophys 393 523 2002 http www mpa garching mpg de Hydro RGRAV index html 153 21 22 23 24 25 26 27 28 29 80 31 82 E E Flanagan S A Hughes Phys Rev D57 4535 4565 1998 N Arnaud et al Phys Rev D59 082002 1999 T Pradier et al Phys Rev D63 042002 2001 N Arnaud et al Phys Rev D67 062004 2003 P Flandrin Temps fr quence Hermes 1993 P Mallat Exploration des signaux en ondelettes Ellispes 1999 N Otsu A Threshold Selection Method from Grey Level Histograms IEEE Trans Systems Man and Cybernetics 9 1 377 1979 J Hoshen R Kopelman Percolation and cluster distribution 1 Cluster multiple labelling technique and critical concentration algorithm PRB14 3438 1976 M Shensa The Discrete Wavelet Transform Wedding the A Trous and Mallat Algorithms IEEE Transactions on Signal Processing 40 10 2464 1992 I Daubechies W Sweldens Factoring wavelet transforms into lifting steps
22. 23 LIVINGSTON expected sensitivity curv parameters static const double S1_LIVINGSTON 1 94e40 static const double S2_LIVINGSTON 1 2e 36 static const double S3_LIVINGSTON 4 68e 46 static const double S4_LIVINGSTON 2 88e 46 static const double FCUT_LIVINGSTON 150 LIVINGSTON geometry angles are in degres static const double LONGITUDE_LIVINGSTON 90 77 static const double LATITUDE_LIVINGSTON 30 56 static const double AZIMUTHAL_LIVINGSTON 333 0 static const double ARM_ANGLE_LIVINGSTON 90 0 GEO 8192 kHz sampling static const int FZERO_GEO 8192 GEO noise SIGMA_NOISE SENS_NOISE sqrt FZERO 2 static const double SIGMA_NOISE_GEO 4 e 21 static const double SENS_NOISE_GEO 4 e 23 GEO expected sensitivity curv parameters static const double S1_GEU 2 18e3 static const double S2_GEO 5 1e 43 static const double S3_GE0 2 0e 46 static const double FCUT_GEO 150 150 GEO geometry angles are in degres static const double LONGITUDE_GEO 9 81 static const double LATITUDE_GEO 52 25 static const double AZIMUTHAL_GEO 158 78 static const double ARM_ANGLE_GEO 94 33 TAMA 8192 kHz sampling static const int FZERO_TAMA 8192 TAMA noise SIGMA_NOISE SENS_NOISE sqrt FZERO 2 static const double SIGMA_NOISE_TAMA 4 e 21 static const double SENS_NOISE_TAMA 4 e 23 TAMA expected sensitivity curv parameters static const double SPEND_
23. Correlators already used in Merlino For the moment the end effect is not treated Hence it is recommended to overlap the data chunks 121 2 11 2 Class Detector The content of the class Detector is given for the sake of completion but one can not instance such an object There is only one public method which returns the name of the algorithm The methods relative to the filter memory facility are protected Include file BuF_Detector hxx Constructors Detector string aName Virtual methods bool Configure string aFileName bool Amplitude Vector lt double gt amp aInput Matrix lt double gt amp aAmplitude int Threshold Matrix lt double gt amp aAmplitude list lt Event gt amp aEventList bool Detect Vector lt double gt S alnput Matrix lt double gt amp aAmplitude list lt Event gt amp aEventList bool End Public method string GetName Non member function int Statistique ofstream amp aFile int aGPS list lt Event gt amp aList The functionalities of the 5 virtual methods may vary according to their specific implementation for each filter but one can extract some common features By now the algorithm configuration parameters are stored in an ASCII file whose name is passed as argument of the the method Configure The syntax of this configuration file depends on the algorithm definition and hence on the specific implementation of the method Configure The Amplitude
24. Frequency cut 5 00 78 AR 512 1 9 0243985366985069e 01 2 5 9914915520694234e 01 511 9 9890829891022292e 04 512 5 2225395422075837e 04 The GetData method fills a vector with Gaussian and colored data followed P yli X angli k Sigma x wli k 1 where w i is a white and Gaussian random variable of variance equal to 1 and az are the AR coefficients The GetPSDFromData method estimates the one sided PSD of a vector of generated data using the AR parameters of the class The GetPSDFromCoefficients method computes the one sided PSD numerically using the AR pa rameters following 1 E ee agen rifk The GetSamplingFrequency returns the Sampling frequency used in the AR process PSD f 2 x Sigma x The GetPolesNumber returns the number of poles number of AR coefficients of the AR process The GetSigma returns the value of the discrete sigma which is used in the AR process The SetSeed allows to set the internal seed number to the given value This method can be called at any instant The SetSeedLocalTime allows to set the seed number used by MT algorithm of a given number which depends on the local time with a precision of 1 s The internal seed number is set to the value given as argument This method can be called at any instant Program exemple include Noise_GaussianAR hxx include DataContainer_Vector hxx include Spectral_Spectral hxx using namespac
25. GetRandomData Vector lt float gt amp aData double GetRandomData double GetRandomData double aXmin double aXmax double GetGaussianData void SetSeed unsigned long aSeed void SetSeedLocalTime unsigned long aSeed void SetDebug bool aValue bool GetDebug This Random_MT class provides facilities to generate a random variable whose values are in between a given segment min 2max The random variable algorithm used in this class is the Mersenne Twister algorithm 12 Its period is 2 x 19937 1 and it is 5 times faster than the previous Random_NR gener ator The internal seed number is initialized to 0 The SetSeed allows to set the internal seed number to the given value This method can be called at any instant The SetSeedLocalTime allows to set the seed number used by MT algorithm of a given number which depends on the local time with a precision of 1 s The internal seed number is set to the value given as argument This method can be called at any instant The GetRandomData returns a random variable whose value lays in the segment 0 1 Note that there is no change of a new sequence seed at each call of the function The GetRandomData aXmin aXmaz returns a random variable whose value lays in aX min aX maz Note that there is no change of a new sequence seed at each call of the function The GetRandomData Vector lt double gt aData fills a vector with random variable whose value lays in th
26. Hanford Hanford LR Go om go Public methods int GetTwoSidedPSD int aSize Vector lt double gt amp aData int GetOneSidedPSD int aSize Vector lt double gt amp aData int GetTwoSidedASD int aSize Vector lt double gt amp aData int GetOneSidedASD int aSize Vector lt double gt amp aData void SetSamplingFrequency int aValue int GetSamplingFrequency double GetLatitude double GetLongitude double Get Azimuthal double GetArmAngle void SetLatitude double aValue void SetLongitude double aValue void SetAzimuthal double aValue void SetArmAngle double aValue void SetDebug bool aValue bool GetDebug The Hanford constructor consider the ITF Hanford_4k detector with its nominal fo 16384 Hz sam pling frequency The three other constructors allow to choose another Hanford detector and or set up a sampling fre quency if needed otherwise the default value is 16384 Hz Note that for the moment the spectrum power is identical for the 3 LIGO detectors The methods are described in the section of the ITF class Program exemple include ITF_Hanford hxx using namespace _BuL int main int arg char argc int i j int size 10000 Vector lt double gt psd_2k size Vector lt double gt psd_4k size Hanford H4k 10000 H4k GetOneSidedASD size psd_4k for i 0 i lt size 1000 i 67 cout cout lt lt lt lt vector ASD lt lt i lt lt lt lt psd_4k i lt l
27. J Fourier Anal Appl 4 3 247 1998 R E Kalman A new approach to linear filtering and prediction problems ASME J Basic Eng 82 35 1960 R G Stockwell L Mansinha R P Lowe Localization of the complex spectrum the S transform TEEE Trans Sig Proc 44 998 1996 154
28. aData double amp aGP S list lt Event gt amp aEventList int End void EnableRefreshPSD bool aValue int UpdateCurrentPSD Vector lt double gt aData int GetNTemplate int GetWindowSize double GetSigmaTemplate int alndex void SetDebug bool aValue bool GetDebug int GetSamplingFrequency int SetSamplingFrequency int aValue int Correlation Vector lt double gt amp aData Matrix lt double gt amp aAmplitude The gaussian templates definition is given in the PC filter datacard The parameters are the minimal and maximal sigma values Omaz Gmin as well as the maximal SNR loss In order to reduce the wrap around pollution when doing correlation using FFT zero padding is used internally Besides to reduce the leakage effect when doing a FFT on a square data window we have used an improved data windowing method Indeed if one simply convolves the data chunk with a moving taper window one finally loose information contained in the data located at the extremity of the window To avoid such an effect we use the same method as the Welch Overlaping Approach usually used for Spectral estimation That is to say the correlation tapered window is moved with a stride equals to half 10The user does not care about adding zero padiing at the end of the input data vector 133 of its size as shown in Figure 2 23 and the PC output amplitude is the combinaison sum in time domain of all the tapered window output
29. aWindowSize int aMovingDelta Vector lt double gt amp aMemory Public methods int SetMemory Vector lt double gt amp aMemory int SetMemory double aValue void SetDebug bool aDebug Value bool GetDebug Virtual method int Compute Vector lt double gt amp aData Vector lt double gt amp aOut The SetMemory methods allow to change the content of the memory of the algorithm at any moment 33 2 4 1 2 Class Kurtosis The Kurtosis estimator measures the degree of peakedness or the degree of fat tails of a distribution of probability density function P x It is defined as pa b2 2 3 Ha where u4 and ua are the fourth and second central moments of the distribution of x which are defined as pa lt a lt a gt gt fe lt zx gt P x dx 2 4 where lt x gt denotes the expectation value The Kurtosis estimator for a Gaussian distribution is 3 Include file Algorithm_Kurtosis hxx Constructors Kurtosis int alnputSize int aWindowSize int aMovingDelta Public methods void SetDebug bool aDebug Value bool GetDebug int Compute Vector lt do4uble gt amp aData Vector lt double gt amp aKurto int Computelterative Vector lt do4uble gt amp aData Vector lt double gt amp aKurto int RecursiveEstimate Vector lt do4uble gt amp aData Vector lt double gt amp aKurto The Compute method compute the Kurtosis quantity using the input Vector and fill the result in the output Vector The Kurtos
30. as with the result of the core collapse simulation performed by Zwerger and Mueller ZM catalogue 17 and Dimmelmeier Font and Mueller DFM catalogue 19 The parameters of each signal type are given when creating the object and can be modified later on if needed For what concerns the normalization of the waveform one has to give a Signal to Noise Ratio SNR value or a distance value but only for core collapse waveforms ZM and DFM When dealing with SNR normalization the type of noise hypothesis is important The choice between white noise or colored noise hypotheses is done at the level of the constructor and can be modified The following table summarizes the different use cases Signal Normalization Noise Constructor Gaussian SNR White Noise I Gaussian SNR Colored Noise IT amp III DampedSine SNR White Noise I DampedSine SNR Colored Noise IT amp Il GaussianSine SNR White Noise I GaussianSine SNR Colored Noise IT amp III ZM SNR White Noise I ZM SNR Colored Noise TI amp V ZM distance White Noise amp Colored Noise II amp IV amp VI DFM SNR White Noise I DFM SNR Colored Noise Ill amp V DFM distance White Noise amp Colored Noise II amp IV amp VI The SNR is defined as 21 pele 1 Ad a haf Ai f Sn 4 Sn where S f is the one sided power spectral density function of the interferometer In the case of a White Noise hypothesis the SNR is simply SNR de aPadf amp hO
31. coefficients according to the format which is described in section 2 7 3 1 The GetSigma alndex method returns the value of the P estimations of the sigma value The Sigma which is used in the AR model is the one corresponding to the index P Note that the alndex value varies from 0 to P The GetAR methods return the values of the coefficients of the AR filter The SetDeltaFrequency and GetDeltaFrequency methods allow to change and to get the value of the frequency step which is used when computing the autocorrelation function The SetPolesNumber and GetPolesNumber methods allow to change and get the number of poles of the AR filter The SetSamplingFrequency and GetSamplingFrequency methods allow to change and get the 83 number of poles of the AR filter The FindOptimalOrder method determine for a give interferomter the optimal number of poles ac cording to the ML criteria to be explained Program example include Noise_AR hxx include DataContainer_Vector hxx include Spectral_Spectral hxx include ITF_Virgo hxx using namespace _BuL int main int arg char argc int poles_number sampling size double f_cut 4 double PSD_scale_factor 1 e40 poles_number 1024 sampling 20000 size sampling 2 AR algo sampling poles_number algo SetDebug false precision de 0 01 Hz Vector lt double gt PSD 100 size on genere la PSD Virgo itf sa
32. corresponding to the 25 DFM signals Al A2 A3 Ad A5 A em 510 1108 5107 1107 BI B2 B3 B4 B5 B 025 05 09 18 40 GI G3 Ga G T 1 325 1 32 1 31 1 30 1 28 Table 2 2 Initial values of the parameters A 8 and T used in the DFM simulation to produce the 25 waveforms of the catalogue The class DFM provides facilities to get the GW amplitude corresponding to the 25 DFM simulated waveforms using equation 2 26 and assuming we are located in th equatorail plane 0 5 The waveforms of the DFM catalogue 20 have been downsampled at 20 kHz Include file BuS_DFM hxx Constructors DFM string aName double aSNR DFM string aName int aDistance DFM string aName double aSNR int aSampling Vector lt double gt amp aPSD DFM string aName int aDistance int aSampling Vector lt double gt amp aPSD DFM string aName double aSNR ITF Detector alTFName DFM string aName int aDistance ITF Detector alTFName Public methods int GetData double aPosition Vector lt double gt amp aData int GetData double aPosition Vector lt float gt amp aData int SetDFM string aName int SetDistance int aDistance 99 The method GetData fills the data vector with the chosen DFM signal which is located at the relative position given by the argument aPosition This position is a fraction of the data vector size If the size of the DFM signal is longer than the vector data size a warning message is printed out and the vec
33. double aA double aX double amp aGln P a x a gt 0 x gt a 1 Note that the erfc function is available in libc math h on Linux RH 6 1 The erfc function is defined as erfela 3 Pat 16 AO T a T a x e tte dt 0 2 et pol Q a x 1 Pla x 8 Tayo The functions Correlation AutoCorrelation and Convolution are identical to the methods of the class Correl excepts that zero padding technics are not available 17 2 2 5 Class Histogram Include file Constructors Destructors Public methods This class provides a very basic implementation of an histogram object At creation one can switch on BuM Histogram hxx Histogram bool aSavelverUnderF low Histogram int Define int aBinNumber double aMinVal double aMaxVal double aWeight int Define double aBinSize double aMinVal double aMaxVal double aWeight double aScale int Add double aVal double GetValue int aBinIndex double GetMean double GetSigma int GetEntries int GetBinNumber double GetBinSize int GetUnderFlowEntries int GetOverFlowEntries double GetWeight double GetScale int WriteFile string aFileName int WriteHistogramFile string aFileName int WriteBinHistogramFile string aFileName void SetDebug bool aValue bool GetDebug the option to save in a vector the value of the underflow and overflow entries The methods Define defi
34. estimation it consists of averaging overlapping data segments WOSA approach Indeed with a overlap greater than 50 the data region which is underweighted at the extremity of one segment is taken into account in the next segment One has to note that the overlapping feature recovers also some information concerning the autocorrelation contained in data values spanning two adjacent non overlapping data segments LI RATS J 09 3 0 8 4 0 7 _ Hann 4 L Haming 0 6 amp Welch 7 Bartlett 0 5 a 04 E 4 03 4 02 O1 de 0 KA L L L 1 L O 0 200 400 600 800 1000 Figure 2 1 Several window functions are usually used For the choice of a particular function there is a trade off between having a principal lob as narrow as possible and the smallest sidelobs amplitude BuL is providing the following tapers which are represented in Figure 2 1 e Hann w 4 1 cos 2 2nt e Haming w 0 54 0 46cos 5 e Welch w 1 1 4 e Bartlett w 1 N 25 2 3 2 Class Taper Include file Spectral_Taper hxx Constructors Taper Taper int aSize Public methods int GetDataHann Vector lt double gt amp aWindow int GetDataHaming Vector lt double gt amp aWindow int GetDataWelch Vector lt double gt amp aWindow int GetDataBartlett Vector lt double gt amp aWindow int GetRawDataHann Vector lt double gt amp aWindow int GetRawDataHaming Vector lt double
35. gt amp aAmplitude double amp aGPsS list lt Event gt amp aEventList int Detect Vector lt double gt amp aData Matrix lt double gt amp aAmplitude list lt Event gt amp aEventList int Detect Vector lt double gt amp aData double zaGPS list lt Event gt amp aEventList int End int SetThreshold Vector lt float gt aThreshold int GetThreshold Vector lt float gt amp aThreshold int GetNWindow int GetWindow int alndex void SetDebug bool avalue bool GetDebug int MeanCompute Vector lt double gt amp aData int aWindow Vector lt double gt amp aOutput The method SetThreshold changes the initial values of the thresholds defined for each window The method GetThreshold returns the value of the thresholds The method GetNWindow returns the number of windows which are used in paralell The method Get Window int alndex returns the size of the window corresponding to the index alndez 131 2 11 7 Class NF The class NF implements the algorithms described in 22 The algorithm is coded in the public method normFilter This function computes the sum of the squared value of N points over a moving window of size n The output result is returned in a vector filled up to bin N n 1 The end of the vector from bin N n up to bin N 1 is filled with zeros The filter output follows a gaussian probability distribution when n is large enough above 30 The filter output is given by i k n 2
36. int int int int int int int int int int int int int int int int int int MedianThresholding Matrix lt double gt amp aMap SigmaThresholding Matrix lt double gt amp aMap const double aTimesSigma OtsuThresholding Matrix lt double gt amp aMap LocalThresholding Matrix lt double gt amp aMap const string aFileName LocalThresholding Matrix lt double gt amp aMap const Vector lt double gt amp aVector LocalThresholding Matrix lt double gt amp aMap PercentileThresholding Matrix lt double gt amp aMap const double aPercent HKclustering Matrix lt double gt amp aMap HKConnectedDimiClustering Matrix lt double gt amp aMap HKNo0WConnectedDimiClustering const Matrix lt double gt amp aMap HKConnectedDimiProfileClustering Matrix lt double gt amp aMap const Vector lt double gt amp aProfil HKNo0WConnectedDimiProfileClustering const Matrix lt double gt amp aMap const Vector lt double HKNo0WConnectedDimiProfileClustering const Matrix lt double gt amp aMap const Vector lt double HKProgressiveClustering Vector lt double gt amp aVect LocalMaxInfo Matrix lt double gt amp aMap LocalMaxInfo Matrix lt double gt amp aMap const int aRange Connexion const double aDilatation ClusterCoincidence void ClusterCoincidence const vector lt Cluster gt amp aCV1 const double aScalei const vector lt Cl void StackClusterEnergy vector lt Cluster gt amp aCV const Matrix lt double gt
37. int size 10000 Vector lt double gt psd_sdt size Vector lt double gt psd_quartz size Vector lt double gt psd_quartz2 2 size 65 Virgo std 10000 std GetOneSidedASD size psd_sdt for i 0 i lt size 1000 i cout lt lt vector ASD lt lt i lt lt lt lt psd_sdt i lt lt endl Virgo quartz ITF Virgo_Quartz 10000 quartz GetOneSidedASD size psd_quartz for i 0 i lt size 1000 i cout lt lt vector ASD lt lt i lt lt lt lt psd_quartz i lt lt endl quartz GetOneSidedPSD size psd_quartz for i 0 i lt size 1000 i cout lt lt vector PSD lt lt i lt lt lt lt psd_quartz i lt lt endl Virgo quartz_bis ITF Virgo_Quartz 20000 quartz_bis GetTwoSidedASD size psd_quartz2 double longitude quartz_bis GetLongitude double latitude quartz_bis GetLatitude double azimuthal quartz_bis GetAzimuthal double arm_angle quartz_bis GetArmAngle cout lt lt Virgo Longitude lt lt longitude lt lt endl cout lt lt Virgo Latitude lt lt latitude lt lt endl cout lt lt Virgo azimuthal angle lt lt azimuthal lt lt endl cout lt lt Virgo arm angle lt lt arm_angle lt lt endl return Q 66 2 6 3 Class Hanford Include file TTF_Hanford hxx int aSamplingFrequency ITF Detector aName ITF Detector aName int aSamplingFrequency Constructors Hanford Hanford
38. is due to the fact that the FFTW plans are determined only once when the object is created The memory allocation of the temporary vector which are needed for each methods is also performed only once in the constructor This feature will speed up the processing time in case of important call number to these spectral transforms Include file Spectral_Spectral hxx Constructors Spectral int aSegmentSize Taper Window aWindow int aSegmentNumber Spectral PSDType aPSDType Destructors Spectral Public methods int Periodogram Vector lt double gt amp aData Vector lt double gt amp aPeriod int PSD Vector lt double gt amp aData Vector lt double gt amp aPSD int aSampling int PSDWOSA Vector lt double gt amp aData Vector lt double gt amp aPSD int aSampling int ASD Vector lt double gt amp aData Vector lt double gt amp aPSD int aSampling int ASDWOSA Vector lt double gt amp aData Vector lt double gt amp aPSD int aSampling int SpectrumPhase Vector lt double gt amp aData Vector lt double gt amp aPhase int SetWindow Taper Window aWindow int aWindowMean int SetSegmentNumber int aSegmentNumber int WriteSpectrumFile Vector lt double gt amp aData string aFileName int aSampling int WriteSpectrumFile Vector lt double gt amp aData string aFileName int aSampling int aFMin int aFMax double FlatenessOfPSD Vector lt double gt amp aPSD double LocalFlatenessOfPSD Vector lt double
39. lt lt endl cout lt lt AIGO Latitude lt lt latitude lt lt endl cout lt lt AIGO azimuthal angle lt lt azimuthal lt lt endl cout lt lt AIGO arm angle lt lt arm_angle lt lt endl return Q 72 2 7 Noise generators of white and colored noise This package includes all classes dedicated for the production of simulated data The noise models which are considered so far are this list will be upgraded e Gaussian and white e Gaussian and colored Two random generators are provided from which are inherated two Gaussian generators Two models of colored noise are provided an AR and ARMA model classes Finally the facilities to produce the AR and ARMA coefficients from a given PSD are provided in this package 2 7 1 Random generators 2 7 1 1 Class Random_NR Include file Noise_Random_NR hxx Constructors Random_NR Public methods void GetRandomData Vector lt double gt amp aData void GetRandomData Vector lt float gt amp aData double GetRandomData double GetRandomData double aXmin double aXmax double GetGaussianData void SetSeed void SetSeedLocalTime void SetDebug const bool aValue bool GetDebug This Random_NR class provides facilities to generate a random variable whose values are in between a given segment min Umar The random variable algorithm used in this class is based on the drand48 function 11 The internal seed number is initialized to
40. of the vector containing the frequency domain wave form is just identical or twice as big as the time domain wave form template size depending on the zero padding option choice The method SetWindowSize allows to re set the value of the vector containing the time domain tem plate wave form This is interesting when the largest template is short is be able to set by hand the vector size of the template in order to have enought frequency resolution Finaly one has to mention that by default the templates will be generated with zero padding If one does not want to zero pad one has to use the TGaussian constructor which allow to set up the private zero padding variable to false Time domain Oo Un TITTY LA potty 500 1000 1500 2000 2500 3000 3500 40 Frequency domain 500 1000 1500 2000 2500 3000 3500 40 o Time domain Zero paddin o Un MA AO A A leds 500 1000 1500 2000 2500 3000 3500 40 Frequency domain Zero padding UU o TT 500 1000 1500 2000 2500 3000 3500 4000 Figure 2 16 Exemple in time domain and in frequency domain of a Gaussian template wave form One can see the different size of the vector depending on the use of zero padding 103 Time domain size N Frequency domain size N See se aan fs 2 ts N Time domain size N Zero padding Frequency domain
41. signal SetSNR 100 signal SetSamplingFrequency 10000 double x 1 4000 95 cout lt lt signal GetValue x lt lt endl 96 2 8 5 Class ZM The Zwerger amp Mueller catalogue gathers a collection of waveforms produced by the Newtonian simulation of the collapse of a rotating stellar core to a neutron star assuming an axisymetry The GW field amplitude emitted during the core collapse is derived from the quadrupole moment By symetry 2D and axisymetry there is only one non vanishing term A A The GW amplitude is then given by 1 1 AR hy z Siro A 2 25 hy 0 The ZM catalogue contain 78 waveforms obtained for different value of the following parameters e A Angular momentum eB f pal ratio between the rotational and the potential energies e G T reduced adiabatic indice The Table 2 8 5 gives the initial values corresponding to the 78 ZM signals AI A2 A Ad AB A cm 510 1108 5107 1107 Bl B2 B3 B4 B5 B 025 04 09 18 40 Gl G2 G3 G4 G5 T 1 325 1 32 1 31 1 30 1 28 Table 2 1 Initial values of the parameters A 5 and T used in the ZM simulation to produce the 78 waveforms of the catalogue The time of bounce is defined as the moment when the central density reaches a first maximum It mainly depends on the initial pressure reduction parameter T and can be described as 4 EC G T iti to where t 20 99 ms and tg 0 63 ms The class ZM provides facilities to g
42. size 2N Zero padding 0 1fs 2 1s 2N Figure 2 17 Representaion of the size of the time domain and frequency domain of the gaussian template Program example include lt math h gt include BuT_TGaussian hxx include lt fstream gt using namespace _BuL int main int narg char argc int i double sigmaMin 1 e 4 double sigmaMax 1 2e 2 double epsilon 1 e 2 int NTemplate bool padding true int size TGaussian TG sigmaMin sigmaMax epsilon padding TG SetDebug true TG SetWindowSize 4096 NTemplate TG GetNTemplate cout lt lt template lt lt NTemplate lt lt endl Vector lt double gt para para TG GetSigmaTemplate for i 0 i lt NTemplate i cout lt lt Template lt lt i lt lt lt lt TG GetSigmaTemplate i lt lt lt lt parali lt lt endl cout lt lt Template size lt lt TG GetWindowSize lt lt endl cout lt lt Generate WF lt lt endl TG GenerateWaveForm 104 cout lt lt Fill vector with WF lt lt endl Matrix lt double gt WF WF TG GetWaveForm return 0 105 2 9 2 Class TDampedSine to be written 106 2 10 Tiler Placement of templates in the parameter space This package provides a tool to place template in a 2D parameter space 2 10 1 Introduction In order to insure a good signal recovery the match filter technique often requires to ru
43. stored in a list 124 Include file Constructors Public methods Non member functions BuF_Event hxx Event int GetTStart int GetTEnd int GetTMax double Get AMax int GetTolerance string GetAlgoName string GetAlgoType int GetTemplateNumero int GetTemplateNumber double GetGPS void SetTStart int aValue void SetTEnd int aValue void SetTMax int aValue void SetAMax double aValue void SetTolerance int aValue void SetAlgoName string aName void SetAlgoType string aName void SetTemplateNumero int aValue void SetTemplateNumber int aValue void SetGPS double aValue void PushBackMicroEvent Event amp aEvent void PushFrontMicroEvent Event amp aEvent void ClearMicroEvent list lt Event gt GetMicroEvent void SetDebug bool aValue bool GetDebug ostream amp operator lt lt ostream amp aStream Event aEvent ostream amp operator lt lt ostream amp aStream Event aEvent void SortStart list lt Event gt amp aList list lt Event gt ClusterBase list lt Event gt amp aList list lt Event gt ClusterCorrelatorTight list lt Event gt amp aList list lt Event gt ClusterCorrelatorLoose list lt Event gt amp aList Some non member functions have been defined The operator lt lt have been defined for an event and a list of event in order to print out in a pretty format an event or a list of events Program example
44. the cluster e Cluster position in dimension 1 estimated with the barycentre see Cluster e Cluster position in dimension 2 e Cluster moment 11 e Cluster spread in dimension 1 e Cluster spread in dimension 2 e Cluster weight against cluster size the weight is summed for a given size to allow mean estimation though without access to sigma e Cluster moment 11 against cluster size similarly to above e Cluster spread in dimension 2 against position in dimension 1 similarly to above It should be noted that this particular output stands for a time spread distribution and is therefore useful only when the map is correctly oriented 144 Some more parameters seem interesting and are to be implemented e Cluster size against dimension 1 and 2 e Cluster energy against dimension 1 and 2 Include file TeF _Statistics hxx Constructors Statistics Statistics const bool aSaveUverUnderFlow const bool aDebug Public methods int int int int int int int int int int int int int SetDimiSet const double aMin const dimblieuBiepIet const Vector lt double gt aV SetDimiWeights const Vector lt double gt aV SetDim2Set const double aMin const double aStep PowerTrigger vector lt Cluster gt amp aClusterList const double aThreshold PowerTrigger list lt Cluster gt amp aClusterMainList vector lt Cluster gt amp aClusterListIn const PowerTriggerNLSpacing list lt Clust
45. to use a mirror condition at both ends of the data block to provide values for all wavelet coefficients instead of a constant maybe zero value There is only one wavelet available at the moment for the a trous transform the cubic b spline whereas Haar and so called 9 7 wavelets are also implemented for the DWT Future additions e Continuous wavelet transform to have a better frequency resolution with an algorithm close to ATrous e Non mirrored ATrous algorithm maybe as an option to ATrousParametrize where the maximum scale is not given by the data length but such as coefficients do not need non provided values before and after the data chunk in time 2 12 3 Class Cluster Include file TeF _Cluster hxx Constructors Cluster Public methods int Empty void int AddPix const Pix amp aPix int AddPix const int aDim1 const int aDim2 const double aVal int SetPixLocalMax const int aNum list lt Pix gt amp GetAnchor void int Absorb list lt Pix gt amp aListe void MoveCluster const int aDimi const int aDim2 void SetCloseToBorder void void SetImageNum const int aNum void SetComposite void void SetPhysical const double aDimiMin const double aDimiStep const double aDim2Min c void SetPhysical const Vector lt double gt aDimiV const double aDim2Min const double aDim2St void SetPhysical const Vector lt double gt aDimiV const Vector lt double gt aDimiW const double a void ScaleEnergy con
46. 2 The WriteOneSidedSpectrumFile and WriteTwoSidedSpectrumFile methods write out in a file the frequency value and the value of a PSD or ASD stored in a vector it corresponds to The GetLatitude method returns the latitude of the Virgo detector The GetLongitude method returns the longitude of the Virgo detector The GetAzimuthal method returns the azimuthal angle of the Virgo detector The GetArmAngle method returns the angle between the two arms of the Virgo detector The SetLatitude method allows to set the latitude of the Virgo detector The SetLongitude method allows to set the longitude of the Virgo detector The SetAzimuthal method allows to set the azimuthal angle of the Virgo detector The SetArmAngle method allows to set the angle between the two arms of the Virgo detector Program example include ITF_ITF hxx using namespace _BuL int main int arg char argc ITF std ITF Virgo double longitude std GetLongitude double latitude std GetLatitude double azimuthal std GetAzimuthal double arm_angle std GetArmAngle cout lt lt Virgo Longitude lt lt longitude lt lt endl cout lt lt Virgo Latitude lt lt latitude lt lt endl cout lt lt Virgo azimuthal angle lt lt azimuthal lt lt endl return 0 64 2 6 2 Class Virgo This class is inherated from the class ITF hence all the public methods of ITF are available in the Virgo
47. 2 21 The full line histogram shows the output of the Mean Filter around a link between 2 data chunks without using the filter memory facility The dashed histogram shows the same output but using the the filter memory facility 120 A simple mecanism to save the filter memory has been implemented in the Detector class It provides an alternative and efficient method to the data chunk overlapping that is necessary to be able to cover all the data stream This filter memory mecanism consists in saving at each call of the Amplitude or Detect method the last part of the data chunk which are needed to be able to compute correctly the filter output of the next data chunk This algorithm is derived from the one implemented in section 2 4 1 The output of a filter at given time stamp t is computed using the input data located between ti y and t where N is the size of the window Hence for the first data chunk the N points are not correct since they are computed with part of 0 One should also mention that this filter memory algorithm takes into account the fact that the filters are implemented with several windows in parallel Finally it must be said that this mecanism imposes to work with constant input data vector size The comparison between the output around 2 data chunks junction of the MF filter using or not the filter memory option is shown in Figure 2 21 Correlator filters PC and DSC We plan to implement the overlap and add nethod for
48. 2 5 1 We give in the following sections few examples of very simple desgined filters using IIR class 2 5 2 1 Low pass filter A low pass filter is realized using one pole with Q 100 and Q 0 in s plane domain The transfer function see Fig 2 4 is then 1 O 1 TITI FERN ica amplitude a eee n PE uc Pa sn ia a ee 4 10 10 E 2 5 0 Edel mal E ie E iia CAD AR EE 2 3 4 1 10 10 10 10 Figure 2 4 Transfer function amplitude and phase of a low pass filter designed using the IIR class The x axis represents frequency in Hz 39 2 5 2 2 High pass filter The transfer function given in Fig 2 5 has been obtained for and IIR filter containing e atriple poles at 2 100 Hz e a triple zero at N 0 Hz The transfer function in s plane domain is jai 1 pun 10 10 1 w 10 ES 10 1 CA 10 a 10 Ss H s _ Y a E amplitude 1 10 10 10 10 E phase 4 E Il a E 1 10 10 10 10 Figure 2 5 Transfer function amplitude and phase of high pass filter designed using the ITR class The x axis represents frequency in Hz 2 5 2 3 Multi poles and zeros filter The transfer function given in Fig 2 6 has been obtained for and IIR filter containing e 3 poles Q Q 10 10 100 10 and 1000 1000 e 1 zero Q Q 300 10 40
49. A doi de lentes AA 115 2 10 7 Class TilerWithEllipses 117 2 10 7 1 A TilerWithEllipses example 118 2 11 BuE burst filters ira a E a Gd DT i gt Ra ne 120 2 11 1 How to treat the end effects between two data chunks 120 2 11 2 Class DEtECtO0 sos s s es 20e eS ee Bo RO as ER ie oe 122 DALI Class Event E SR Ee A i a ee Be OEE ee en 124 2 13 2 14 2 15 2 16 2 17 2 11 4 Class Example or how to code a new Burst algorithm 128 2 11 05 Class AMES a eee See aes wh oh eA and A GEM ti BUR EIS ye 130 QA T iG Class ME hex A ete i Gh oS aes sgt Nae A aii Beste 131 2117 Class NP ois tye is ae BR RS Seek Be ee EE A 132 2 ES CASS O he sei cette tees Aste ate Se vases Gis Re ante poe poy Cat es oe ge ME EE AERO OO A eNO 133 2 019 Class DS Gris Sy By Syste MR a Gow Bee Faw ee Sk ded Deere 135 II Clas LEY dag doh tech oth be eho A Ad E ote te Oo eta thes age 135 QA TT lassGDE 23 ie Se Se Se Soh ae Nee ee ee we re RRR eB een hc Ke Se Re ete BEE 135 EDGR gece Ree pad eal ae EOE Ee ce ew Se ewe A ORS NE ee ES a aaae a N 136 DD O ONS 225 cn eee A Shen be Gee Awe D ni a ig 136 2 12 9 Class Transforms ss ca aoe a tes Sas de WO nal a Ge wt edge do Ge Hee a Se 136 2 12 3 Class Cluster sd se se a ea EE vin WR EOE ee a A eee Se 52008 139 2 12 4 Class Image ee Eee LR A fe ere ee dre Bea NE See Be 141 2 112 51 Class Statistics staat octets as ees Ee
50. BuL Virgo Gravitational Wave Searches Library for Burst Sources User s Manual v3r4 7th January 2005 VIR MAN LAL 7400 100 N Arnaud M A Bizouard F Cavalier A C Clapson M Davier G M Guidi P Hello N Leroy T Pradier A Vicer Laboratoire de l Acc l rateur Lin aire INFN amp Universita di Urbino Contact M A Bizouard mabizoua lal in2p3 fr Contents 1 Introduction 2 Teds Generalities A SUA Eyes Geechee ashes LSE o a a Bele ge OI 12 Bulb library installation s x04 ay GP a aes Sie eet de RE A Ww ww Ae ws 13 CF eode stylet ack eet A hoe Se Bd dct it M AN MERDE tee EX Namespace irin u sa arr Stee ee ee pd eee ee a e Bul packages description 2 1 Name convention i335 alas foe eh A Soe A oS eee bo da fe 2 2 BuM mathematics for bursts search A AG NO E RR NN 22d Fourier transformer GA a ROR dure ee Oe ee 22 1 2 EFTW description y soo cot a 3 a eon ue AT te bete a kG ee 50 2 2 1 3 Class FFT description 2 22 Glass Corre na rar ee gee Gee Ph A ek RP o ek de Eee tt 2 2 3 Ulass GPS Times 5 2506 hee we EEE ee ee eS OR 2 2 4 Mathematical functions 2 2 5 Glass HISTOST AM a a an at er Daag ak MLM Pat ae 2 2 6 Class CumulativeHistogram 227 O E A Bo ed Ra ete ty en Soe Se ees E Aa AE ARA EA A E Aa EN Deh SOP OCEBAM se a Best teal A ak Ets alec RA ae eee By
51. ComputeSpaceParamArea double aAccuracy void SetStartingCenter double aCoord TilerPoint amp GetStartingCenter void CentersToFile std string aFileName The constructor Tiler int aDim 2 creates a Tiler object for a parameter space of dimension aDim By default it works in a 2D space The virtual Tile method is the one which will do the tiling for each derivated classes The GetDim method returns the dimension of the parameter space For debug prints the SetPrintPeriod int aPeriod method sets the period of print in the choose of centers The program will print the following lines each time aPeriod centers have been tested in the tiling procedure 6000 Centers chosen Total Time Elapsed 137 Since Previous 14 115 It also gives the time elapsed since the start of the program and the time since the previous print By default the PrintPeriod is set to 1000 The GetPrintPeriod method returns the Print Period defined above The SetDecreaseSurfaceFactor double aValue method set the DecreaseSurfaceFactor variable Each time the remaining surface to be tiled divided by the total parameter space surface is decreased by this factor the following print is done Remaining Surface 2 07861 Used Ellipses 359 Total Time Elapsed 18 Since Previous 6 The print gives the remaining surface to tile and the number of ellipses already placed It also shows the time elapsed since the start of the program and the t
52. DWOSA has been already performed Otherwise it returns false 30 Program example include Spectral_Spectral hxx using namespace _BuL int main int narg char argc int i double x double sigma 10 int segment_size 200000 int segment_number 10 double khi int size segment_size segment_number Spectral titi segment_size Taper Hann segment_number Spectral TwoSided Vector lt double gt aaa size Vector lt double gt bbb segment_size for i 0 i lt size i x i size 8 5 aaali exp x x 2 sigma sigma titi PSD aaa bbb 20000 khi titi FlatnessOfPSD bbb cout lt lt flatness of PSD lt lt khi lt lt endl titi WriteSpectrumFile bbb psd dat 20000 titi PSDWOSA aaa bbb 20000 khi titi FlatnessOfPSD bbb cout lt lt flatness of PSD WOSA lt lt khi lt lt endl titi WriteSpectrumFile bbb wosa dat 20000 return 0 31 2 4 Algorithm algorithms for statistical analysis 2 4 1 Slicing window algorithms A lot of algorithms are estimating a quantity over a chunk of data using a moving window by step of one or several bins Such algorithms need to keep in memory the value of the last data points in order to avoid discontinuities in the algorithm response at each new data chunk The class MovingSlice has been designed to implement these facilities which are common to all the slicing window algorithms The o
53. Frequency step lt filterEndFrequency step filterStep cout lt lt tested frequency lt lt step lt lt endl filtre Setup samplingFrequency step loopSize filtre Filter data line window filtre Filter data line for k 0 k lt dataSize loopSize 1 k result k index line k resu_1 index line 0 resu_2 index line 10 resu_3 index line 20 cout lt lt index lt lt endl index data WriteFile testdataGoert dat line WriteFile testlineGoert dat result WriteFile testsweepGoertzel dat resu_1 WriteFile toto dat resu_2 WriteFile titi dat resu_3 WriteFile tutu dat return 0 61 2 6 ITF detector characteristics This package is intended to contain all the header files defining parameters which are specific to the detectors For the moment it includes all the Earth interferometers header files e GEO ITF_GEO_Constant hxx e Hanford ITF_Hanford Constant hxx e Livingston ITF Livingston_Constant hxx e TAMA ITF_TAMA_Constant hxx e Virgo ITF_Virgo_Constant hxx e AIGO ITF_AIGO_Constant hxx Their content is given in 2 16 A generic class called ITF is inheritated by several classes which correspond to each Earth based inter ferometers in activity The class ITF provides methods to get mainly geometrical information about an interferometer while the inhnerated classes provide a mean to get th
54. In that way all points have identical weight in the amplitude output Finalyy by default the PSD is estimated using the input data each time the amplitude is computed This default behavior can be changed at any moment for instance when the input data are not good enough as just after a de lock period In that case the PC algorithm will used the previously estimated PSD Tapered data window Hann Data chunk overlap region Input data Output data ET ER A dame Figure 2 23 The Welch Overlapping Approach is used to compute the correlation between small data window and Gaussian template This method has the advantage to give an identical weight to all the input data points as shown in lower part of the figure An example of a PC datacard is given below Example PARAMETERS DEFINITION ALGO Peak Correlator CLUSTER CorrelatorLoose FILTERMEMORY No SIGMA 1 0e 4 1 2e 2 en seconde EPSILON 1 0e 2 WINDOW 16384 THRESHOLD 5 6 6 The Configure method computes the number of templates and generates the wave forms in the 134 Fourier space which are stored in memory The Amplitude method calculates for all the templates the correlation function between the templates and the data vector using the Fourier Transform of both The filter output is given for template i using the Fourier Transfrom by 1 Fate F gausso VTfs0i PSD f Note that the filter output is normalized in ord
55. N 2 H s gt x gt 2 11 o 8 si Ms si 8 0 Grouping poles in complex conjugate pair the product factor if eq 2 10 writes 1 gt HA 5 2 12 82 20 Re sx 02 2 12 Once we have determined the poles the low pass Butterworth filter can be implemented using the IIR class cascading second order ITR filters defined by Q Q Using Eq 2 12 it is easy to see that with reference of convention defined in section 2 5 1 one has O Qe 2 13 Qn Seen 2 14 In case N is odd the single pole factor gives w Q Q 0 2 16 E Rh Besides a high pass filter can be constructed from a low pass filter by means of a transformation Let s consider the low pass transfer function _a0 aZ a Hp AAA A 2 17 Po yd bye 10 To get the corresponding high pass filter let 1 cos wipA twnpA pea o o 2 18 1 02 cos 222P2 3 hp Expanding Eq 2 17 one gets for the high pass filter coefficients ag 20 09 taza 6 1 bia boa a aces al atia b 2a bi 1 a7 2b20 a aoa ayataa bo batb Choosing a 0 implies Whp wipA T w If 2 19 and ao ao a a by b1 ag a2 by da Butterworth band pass filters are then built by cascading a high pass and a low pass filter on data We give in the following section some example of low high and pass band filters and their effect on a colored noise sequence 47 2 5 4 1 Butterworth low pass filter Low pass band frequency 100 Hz Order N 20
56. Pal 2 o2 On 00 n 00 89 2 8 1 Class Source The class Source is a generic class containing virtual methods which are implemented in the derived classes and hence one cannot instance a Source object Nevertheless the Source class contains some public methods which are directly inherated by the derived classes This is the reason why we give the list of the public methods of the Source class in this user s manual When one creates a specific source object one has to specify the noise model which is supposed to be used using the appropriate constructor Default noise caracteristics when a Source object is created e White Noise The sampling frequency and the white noise sigma are extracted from the Virgo_Constant hxx header file which values are given in sections 2 7 2 and 2 6 e Colored Noise One gives the PSD in an input Vector or one gives the name of the interferometer whose PSD and sampling frequency are taken from 2 16 Here is a list of the methods which can be used in any specific source class Public methods int SetDebug bool aValue bool GetDebug int GetSamplingFrequency int SetSamplingFrequency int aValue int SetSigmaNoise double aSigma double GetSigmaNoise int SetSNR double aSNR int SetPSD int aSampling Vector lt double gt aPSD int SetITF ITF Detector aName e The following methods have action only in the case of a white noise The method SetSamplingFrequency changes t
57. R int aSampling int aPolesNumber Public methods int ComputeAR Vector lt double gt amp aPSD double aPSDScaleFactor double aMinimalFrequency int ComputeAR ITF Detector aITF double aMinimalFrequency int SaveAR string aFileName double GetSigma int alndex double GetAR int alndex void GetAR Vector lt double gt amp aAR int SetDeltaFrequency double aValue double GetDeltaFrequency int SetPolesNumber int aValue int GetPolesNumber int GetSamplingFrequency int SetSamplingFrequency int aSampling int FindOptimalOrderMDL ITF Detector aITF double aMinimalFrequency int aDataSize void SetDebug bool aValue bool GetDebug The ComputeAR method compute the AR coefficients from the PSD given as argument using the Durbin Levinson regression to solve the Yule Walker equations 13 The aPSDScaleFactor argument is the scale factor which has been applied to the PSD in order to avoid to manipulate very small number The aMinimalFrequency argument is the frequency below which the PSD is taken as constant and equal to the value at aMinimalFrequency The ComputeAR method compute the AR coefficients taking into account the official PSD which corresponds to the given interferometer The aMinimalFrequency argument is the frequency below which the PSD is taken as constant and equal to the value at aMinimalFrequency The SaveAR method saves in the file whose name is given as argument the AR filter
58. STime constructor To summarize the GPS time conversion is since January 1 1999 given by GPS UTC GPS_EPOCH 13 seconds Include file BuM_GPSTime Constructor GPSTime Public methods timet GPSToSys time_t aGPSTime time_t SysToGPS time t aSysTime timet GPSToUTC timet aGPSTime timet UTCToGPS timet aUTCTime string TimeToString time t aTime string GPSToUTCTimeString timet aGPSTime string GPSToSysTimeString timet aGPSTime string DateToString timet aUTCTime string GPSToUTCDateString timet aGPSTime string GPSToSysDateString timet aGPSTime void SetDebug bool aDebug Value bool GetDebug 15 2 2 4 Mathematical functions Include file BuM Eunction hxx Functions void Integral Vector lt double gt aData double aStep double amp aSum void Integral Vector lt float gt aData float aStep float amp aSum void IntegralNorm Vector lt double gt aData double aStep double amp aSum void IntegralNorm Vector lt float gt aData float aStep float amp aSum void IntegralModule Vector lt double gt aData double aStep double amp aSum void IntegralModule Vector lt float gt aData float aStep float amp aSum double GetThresholdKhi2 double aFalseAlarm int aN double GetThresholdGaussian double aFalseAlarm double aSigma int Correlation Vector lt double gt amp aDataOne Vector lt double gt amp aDataTwo Vector lt double gt amp aData int Correlation Vector lt complex lt double gt gt amp aDataOne
59. Se tes A D 2 3 1 Non parametric spectral estimation Dard Glass Taper ti toi Utes ete he oe ela an cas tee he E E D Be At ds ae 2337 Class Spectral iiss sees ee ee DOR ed PE Re a Se ge te 2 4 Algorithm algorithms for statistical analysis 2 4 1 Slicing window algorithms 2 4 1 1 Class MovingSlice 2412 Class Kurtosis es cae AR PLU EM NS ones a ee ee A 2 4 2 Mohanty s algorithm 2 5 LinearFilter time domain linear filters 2 01 Introductions o ee eee Se ae Te ee AS 2 5 2 How to design a linear filter with class IIR using the zeros and poles design 2 5 3 2 5 4 2 5 5 2 5 6 25 21 bow pass Titer sio Bk a SOP eus Bia eee oe es 2 02 32 High pass filter ss eae se de Ga A eR Roe Ee A 4 SOS 2 5 2 3 Multi poles and zeros filter Class MR at tia dhe E pto Minne Stony ol ha IIR filter using the Butterworth approximation 2 5 4 1 Butterworth low pass filter 2 5 4 2 Butterworth high pass filter 2 5 4 3 Butterworth pass band filter IIR filter using the Chebyshev approximation Class B tteworth 4 dai a A Re Oe A eee 20 7 Class Matchhilter s grele ma fe uh A ee eee a ete
60. TAMA 7 68e 32 static const double SMIRROR_TAMA 5 2e 43 static const double SSHOT_TAMA 9 e 46 static const double FCUT_TAMA 400 TAMA geometry angles are in degres static const double LONGITUDE_TAMA 139 54 static const double LATITUDE_TAMA 35 68 static const double AZIMUTHAL_TAMA 315 0 static const double ARM_ANGLE_TAMA 90 0 AIGO AIGO geometry angles are in degres static const double LONGITUDE_AIGO 115 70 static const double LATITUDE_AIGO 31 35 static const double AZIMUTHAL_AIGO 90 unknown value static const double ARM_ANGLE_AIGO 90 0 151 2 17 Appendix 3 alblgl a2blgl a3blgl adblgl alblg2 a2blg2 a3blg2 adblg2 alblg3 a2blg3 a3blg3 a4blg3 alblg4 a2blg4 a3blg4 adblg4 alblg5 a2blg5 a3blg5 a4blg5 alb2g1 a2b2g1 a3b2gl a4b2gl alb2g2 a2b2g2 a3b2g2 adb2g2 alb2g3 a2b2g3 a3b2g3 a4b2g3 alb2g4 a2b2g4 a3b2g4 a4b2g4 alb2g5 a2b2g5 a3b2g5 a4b2g5 alb3gl a2b3g1 a3b3gl adb3gl alb3g2 a2b3g2 a3b3g2 a4b3g2 alb3g3 a2b3g3 a3b3g3 a4b3g3 alb3g4 a2b3g4 a3b3g4 a4b3g4 alb3g5 a2b3g5 a3b3g5 a4b3g5 a2b4g2 a3b4g2 a4b4g2 a2b4g3 a3b4g3 a4b4g3 azb4g4 a3b4g4 a4b4g4 a2b4g5 a3b4g5 a4b4g5 a2b5g4 a3b5g4 adbdg4 a2b5g5 a3b5g5 adb5g5 Table 2 3 Table of the names of the 78 ZM waveforms alblgl a3blgl a4blgl a4b1g2 alb2g1 a3b2g1 a3b2g2 adb2g2 a4b2g3 a3b2g4 alb3gl a3b3gl alb3g2 a3b3g2 alb3g3 a3b3g3 alb3g5 a3b3g5 a2b4gl a3b4g2 a4b4g4 adb4g5 a3b5g4 a4b5g4 a4b5g5 Table 2 4 Table of the names o
61. Tiler Point amp aPoint method allows to define the function which returns the value of the XYCoeff for the ellipse centered on aPoint The SetSafetyFactor double aValue method sets the safety factor to aValue The GetSafetyFactor method returns the value of the safety factor The SetAngleOffset double aValue method sets the angle factor to aValue The GetAngleOffset method returns the value of the angle offset The SetTurnDirection double aValue method sets the turn direction If a Value is greater than zero TurnDirection is anticlockwise and clockwise if not The GetTurnDirection method returns the turn direction 1 means anticlockwise and 1 clockwise 117 The SetMinimalCoveredSurface double aValue method sets the minimal covered surface criteria to aValue The GetMinimalCoveredSurface method returns the value of the minimal covered surface criteria The SetNbLoopsComputeCoveredSurface int aValue method sets the number loops used to compute the surface covered by an ellipse in the parameter space to aValue The GetNbLoopsComputeCoveredSurface method returns the previously defined variable The Tile method performs the tiling of the parameter space It returns the number of used ellipses The TotalEllipsesArea method returns the sum of all ellipse areas It should be call after the tiling The EstimatedCoveredArea method returns the parameter space area covered by the ellipses The WriteKu
62. ace _BuL int main int narg char argc int i double temps int taille 1000 double frequency 3000 0 double sigma 2 0 Vector lt double gt signal taille 0 0 Histogram gaussienne false gaussienne SetDebug true for i 0 i lt taille i temps double i taille 2 double frequency signal i 100 0 sin 2 M_PI i taille 2 frequency signal WriteFile titi dat gaussienne Define 100 0 100 1 0 cout lt lt taille histo lt lt gaussienne GetBinNumber lt lt endl cout lt lt remplissage histo lt lt endl for i 0 i lt taille i gaussienne Add signal i cout lt lt valeur a indice 30 lt lt gaussienne GetValue 30 lt lt endl cout lt lt valeur a indice 0 lt lt gaussienne GetValue 0 lt lt endl cout lt lt valeur a indice 100 lt lt gaussienne GetValue 100 lt lt endl 19 gaussienne WriteHistogramFile toto dat gaussienne WriteFile tutu cout lt lt mean lt lt gaussienne GetMean lt lt endl cout lt lt sigma lt lt gaussienne GetSigma lt lt endl return 0 20 2 2 6 Class CumulativeHistogram Include file BuM_CumulativeHistogram hxx Constructors CumulativeHistogram Destructors CumulativeHistogram Public methods int Define const int aBinNumber const double aMinVal const double aMaxVal int WriteSigmaHistogramFile
63. after normalizing both mean and standard deviation using the call values Note that is normalization can be done using Import e MedianThresholding cuts at the pixel value median that is mean between minimum and maximum values e SigmaThresholding cuts at the given fraction of the estimated sigma of the pixel values o lt E lt E gt gt e OtsuThresholding uses the Otsu algorithm 27 to define the threshold This relies on the hypothesis of a bimodal distribution of the pixel energy to isolate the higher values part e LocalThresholding methods are not fully supported and try to estimate some zonal threshold 143 e Percentile Thresholding uses the pixel energy histogram that must be constructed beforehand with BuildEnergy Histo It then sums the histogram values till reaching the asked fraction of total pixel number This method currently is meaningless on non positive maps see BuildEnergyHisto HK clustering 28 is in charge of the identification of pixel clusters after determination of the pixel value threshold It uses a strict definition of a cluster a group of connected pixels with a 8 connectivity rule horizontal vertical and diagonal neighbors It overwrites the map with zero for pixels below the thresh old and a cluster number for those above Its main interest is filling a list of Cluster objects thus achieving a first data reduction step Several implementations are proposed They achieve cluster connection ov
64. age of power from one portion of the PSD to another Indeed the sidelobs convolution involves frequencies which are distant from the one of S f The estimated PSD is biased at frequencies for which S f is small compared to the rest of the PSD When N gt x F f f f f lt Sy f gt is hence an aymptotically unbiased estimator of S f F f NDK E In practise one needs a mean to reduce the sidelobs effects which are especially critical for large dynamic range spectrum estimation One technique tapering consists in multiplying in time domain bin by bin the discrete data set x by a window w function data taper and estimate the PSD of this new data set fo 2 lt Sn f taper H t f S f af fo 2 where e H f 5 wyermift 2 t 1 The key idea behind tapering technique is to select the window function w such that H f has much smaller sidelobs than F f This leads to select taper which regularize the window at the extrema see Figure 2 1 Nevertheless this regularization induce an increase of the estimation variance due to the loss of data information at the extrema of the data analysis window This drawback can be overcomed by performing an average estimation of the PSD Sthe other is pre whitening which reduce the dynamic range 24 The most simple is to average the power spectrum using several data blocks or the same but with a reducing frequency resolution The Welch approach 7 gives a better
65. amp aMap const int alnt int int CoherenceHK const Matrix lt double gt amp aMapA const Matrix lt double gt amp aMapB1 const Matrix lt CoherenceHK const Matrix lt double gt amp aMapA const vector lt Matrix lt double gt gt amp aMapB1 cons void Set0ffsetBAdim2 const int aVal int int int int int int int int int int GetClusterNb void GetPixOnNb void GetLocalMask Matrix lt double gt amp aMap GetClusterList vector lt Cluster gt amp aClusterList GetLastClusterList vector lt Cluster gt amp aClusterList FlushMapClusters void FlushMapConnections void WriteEnergyHisto const string aFileName WriteMapFile Matrix lt double gt amp aMap const int aMaxDimi const int aMaxDim const str WriteClustersFile const string aFileName void SetDebug const bool aValue bool GetDebug 142 Some methods are not fully supported yet and therefore best left alone These are Connexion Cluster Analysis and Tests This class is supposed to provide image analysis tools with a focus on time frequency needs It imports TF maps allegedly built with Transforms to extract information The end product is very restricted at the moment a n tuple file of cluster parameters a fixed set of histograms Future development should probably give some leeway on the class output and maybe offer some trigger method unless a dedicated Trigger class is deemed necessary Nota bene some analysis tools
66. are sensitive to the axis signification provided through the integer aTrans form Besides other methods have only been implemented for a given orientation indicated in their name At construction the dimension of the map a Matrix object is asked though it could be recovered Internal objects are allocated memory in relation to these dimensions which can be changed without having to create a new object via Resize Map import can be done from a Matrix with Import or an ASCII file using ReadMapFile It should be noted that the map is overwritten by several algorithms most notably the clustering one Setting the aDebug parameter to true makes most methods overwrite the map The first possible output is the histogram of the map pixel energy built with BuildEnergyHisto and saved to a file by WriteEnergyHisto Note that some transforms are non positive thus making this histogram unreliable at the moment since separation of negative and positive values is not implemented Are also available the mean and standard deviation profiles computed with respect to time as defined by the map orientation Map analysis starts with thresholding Several methods have been implemented without knowledge of their adequacy to our data type e Thresholding does hard thresholding at the given value previous knowledge of map values is rec ommended e MeanThresholding cuts at the mean pixel value e MeanAndSigmaThresholding cuts at the given value
67. ate objects e The BuL library avoids as much as possible the copy of objects Most of the time the classes methods API make use of the objects address 1 4 Namespace The BuL library defines a namespace _BuL which embed all BuL classes and functions It is then mandatory to declare the use of the BuL namespace in any user s program Chapter 2 Bul packages description In the following the different classes and functions are described The documentation is reduced to the public users methods in order to give a rapid overview of the present functionalities For a complete ref erence especially for the private fields and exact types reference one will refer to the reference document provided by Doxygen 2 On the contrary the necessary include file names are systematically given 2 1 Name convention e fo sampling frequency e S f One sided Power Spectral Density function also named PSD The square root of S f is the One sided Amplitude Spectral Density ASD 2 2 BuM mathematics for bursts search 2 2 1 Class FFT 2 2 1 1 Fourier transform The class FFT has been defined in order to provide a friendly interface to a lower level FFT algorithm library This interface allows use of the complex container defined in the C DCPP package This class is an interface to FFTW library We have chosen the following definition of the forward Fourier Transform FT of a function h t h f f sat The FFT class allows to perfo
68. bject is created The memory allocation of the temporary vector which are needed for each methods is also performed only once in the constructor This feature will speed up the processing time in case of important call number to these spectral transforms Include file BuM_Correl hxx Constructors Correl int aSize bool aZeroPadding Public methods int Correlation Vector lt double gt amp aDataOne Vector lt double gt amp aDataTwo Vector lt double gt amp aData int Correlation Vector lt complex lt double gt gt amp aDataOne Vector lt complex lt double gt gt amp aDataTwo Vector lt complex lt double gt gt amp aData int AutoCorrelation Vector lt double gt amp aDataOne Vector lt double gt amp aData int AutoCorrelation Vector lt complex lt double gt gt amp aDataOne Vector lt complex lt double gt gt amp aData int Convolution Vector lt double gt amp aDataOne Vector lt double gt amp aDataTwo Vector lt double gt amp aData int Convolution Vector lt complex lt double gt gt amp aDataOne Vector lt complex lt double gt gt amp aDataT wo Vector lt complex lt double gt gt amp aData void SetDebug bool aValue bool GetDebug One can apply the different methods on zero padded data vectors or on normal data vectors In the case of zero padding the dimension aSize which appears in the Correl constructor is the physical size of the data vector and hence the input vectors size must be eq
69. ble gt data dataSize Vector lt double gt line dataSize int k double lineFrequency1 100 0 double lineFrequency2 101 0 double filterStartFrequency 90 0 double filterEndFrequency 110 0 double filterStep 0 01 int loopSize 20000 int samplesize 20000 double amplitude 10 double step int index double temp int loop_nb filterEndFrequency filterStartFrequency filterStep cout lt lt data size lt lt dataSize lt lt endl cout lt lt loop_nb lt lt loop_nb lt lt endl Vector lt double gt window samplesize Vector lt double gt resu_1 loop_nb Vector lt double gt resu_2 loop_nb Vector lt double gt resu_3 loop_nb Matrix lt double gt result int dataSize loopSize 1 0 int filterEndFrequency filterStartFrequency filterStep 1 0 0 0 noise generators 60 Gaussian_MT NoiseGen 1 0 NoiseGen GetData data add lines for k 0 k lt dataSize k temp amplitude sin 2 0 M_PI k lineFrequencyi samplingFrequency temp amplitude sin 2 0 M_PI k lineFrequency2 samplingFrequency M_PI 3 0 data k temp filtre Setup samplingFrequency lineFrequency1 loopSize cout lt lt one value lt lt filtre Filter data lt lt endl Taper hann samplesize hann GetDataHann window test windowing effect data 1 0 sweep filter frequency index 0 for step filterStart
70. blic methods int Forward Vector lt complex lt double gt gt amp input Vector lt complex lt double gt gt amp output int Backward Vector lt complex lt double gt gt amp input Vector lt complex lt double gt gt amp output int Forward Vector lt double gt amp input Vector lt complex lt double gt gt amp output int Backward Vector lt complex lt double gt gt amp input Vector lt double gt amp output int Forward Vector lt double gt amp input Vector lt double gt amp output int Backward Vector lt double gt amp input Vector lt double gt amp output void SetCheckVectorSize bool aValue int ExportLastWisdomFile char fileName int SaveWisdomFile char fileName int GetSize void SetDebug bool debugValue bool GetDebug int WriteFile Vector lt complex lt double gt gt amp aData string aFileName int aSampling FFT Output aPart int WriteFile Vector lt double gt amp aData string aFileName int aSampling int WriteFile Vector lt complex lt double gt gt amp aData string aFileName int aSampling FFT Output aPart double aFMin double aFMax int WriteFile Vector lt double gt amp aData string aFileName int aSampling double aFMin double aFMax In order to set up the FFTW options one has defined several enumeration types FFT Type FFT Wisdom FFT Output and FFT Direction Those enumeration types are defined inside the class FFT We give their definition below and in sec
71. coefficients which correspond to the PSD of the data he wants to simulate to know how to determine such coefficients see section 2 7 5 In the second constructor II the user gives the name of the interferometer wich corresponds to the data he wants to generate All the ARMA coefficients files corresponding to each official interferometer PSD curves have been produced The format of the ARMA coefficients file is the following A comment line is indicated with Info The Sigma value is the discrete sigma value The Sampling value gives the sampling frequency of the process The Frequency info line gives the minimal frequency which is simulated Below this frequency the spec trum is flat The AR value gives the number P of poles of the ARMA filter The AR coefficient values are then given The MA value gives the number Q of zeros of the ARMA filter The MA coefficient values are then given Sigma 1 8522732942547901e 20 Sampling 20000 Info Frequency cut 10 00 AR 38 1 3 3181070132757395e 01 2 3 7689826182215138e 00 38 2 6785776758171774e 11 MA 16 0 8 5722580329088460e 01 1 1 0604343586542779e 00 81 2 2 7890994776648230e 00 16 2 3504054468144381e 02 The GetData method fills a vector with Gaussian and colored data followed Q yli X axyli k Sigma x S be wf i k 1 k 0 where wii is a white and Gaussian random variable of variance equal to 1 az are the AR coefficients and by are t
72. const int aMother ATrous Vector lt double gt amp aData Matrix lt double gt amp aMap DyadicParametrize const int aMother Dyadic Vector lt double gt amp aData Matrix lt double gt amp aMap void BuildFrequencyProfile Matrix lt double gt amp aMap Vector lt double gt amp aProfile void BuildFrequencyProfileAndSigma Matrix lt double gt amp aMap Vector lt double gt amp aProfile Vec int int int GetDimi void GetDim2 void GetOverlap void double GetAbsMax void int int GetDimi0ffset void GetDim20ffset void double GetDimiMin void double GetDimiStep void double GetDim2Min void double GetDim2Step void int int int SetStartTime const double aTime GetTimeLength GetFrequencyLength double GetTimeMin void double GetTimeStep void double GetFrequencyMin void double GetFrequencyStep void void GetMapTimes Vector lt double gt amp aVector void GetMapFrequencies Vector lt double gt amp aVector void GetMapSigmas Vector lt double gt amp aVector double GetNLSpacingCoef void void WriteTransformWindow const string aFileName void WriteMapTimes const string aFileName void WriteMapFrequencies const string aFileName int WriteMapFile const Matrix lt double gt amp aMap const int aMaxDimi const int aMaxDim con void SetDebug const bool aValue 137 bool GetDebug Each transform method comes with a parameterizing method that sets up the r
73. culate the amplitude of the filters to apply the thresholds on the amplitude and finally to cluster the different triggers output to define the events All the Burst algorithms classes which are described im the next sections derives from the class Detector and implements the virtual methods of the class Detector Each class algorithms is provided with a test application In addition to this it contains a class DataCard whose role is to provide methods to read easily a data card Besides we have given a simple definition of an event in the search of Burst context which leads to the definition of a class called Event 2 11 1 How to treat the end effects between two data chunks Since we have to break the data stream up into data chunks the filter output will not be correct if one does not care about the fact that the beginning of each output is influenced by the last input data of the previous data chunks Two technics are implemented to cope with end effects depending on the nature of the filters Finally it should be mention that the Burst filter classes offer the possibility to enable and disable the specific data chunk end treatment This is defined at the Configuration phase e Time domain filters ALF NF and MF end of data chunk input size window size window size Ne ET 2 filter memory O filter memory OFF 3 pe ils 24 pde ol Ll Ll i 200 400 600 800 1000 Figure
74. d of filter is always stable e TIR Infinite Impulse Response the impulse response has an infinite number of non zero values An IIR filter is stable if h n a u n with la lt 1 and ufn is the unit step discrete function It is usually easier to figure out the filter response in terms of zeros and poles representation instead of using the impulse response We consider in the following the sub class of LTI such as 5 azy n k 5 by a n k k 0 k 0 We will consider without loss of generality that ao equals to 1 Such a system is not absolutly causal linear nor time invariant Using the z transfrom the difference equation becomes N M 5 akz Y z y biz X 2 k 0 k 0 The z transfrom of a sequence x n is the analog in discrete time domain of the Laplace transform 2 transfrom CO X z y a njz k Laplace transform s testa The tranfer function of the system in the z plane is Y 2 _ Veto biz X z ao agz H z The transfer function can be factorized and one obtains a fractional equation bo eci 1 x27 ao Tli 1 1 dea Tit should also be mentioned that the method of filter design by impulse invariance suffers from aliasing 8bilateral definition 35 where cy are the M complex zeros and d are the N complex poles of the tranfer function Note that the FIR case corresponds to a filter where all a are null except ao H z has no pole except in z 0 An IIR filt
75. dat psd ASD data out int sampling psd WriteSpectrumFile out spec dat sampling filter Filter data output output WriteFile filtered_data dat psd ASD output out int sampling psd WriteSpectrumFile out spec_fil dat sampling gene GetData data data AppendFile data dat filter Filter data output output AppendFile filtered_data dat return 0 55 2 5 7 Class MatchFilter This class provide methods for performing the match filtering between a data vector and a given template 56 Include file LinearFilter MatchFilter hxx Constructors MatchFilter int aSampling Vector lt double gt amp aPSDTowSided Vector lt double gt amp aTemplate Destructors MatchFilter Public methods int SetPSD Vector lt double gt amp aPSDTowSided int SetTemplate Vector lt double gt amp aTemplate int SetTemplateAndPSD Vector lt double gt amp aTemplate Vector lt double gt amp aPSDTowSided int Filter Vector lt double gt amp alnput Vector lt double gt amp aDutput double GetNormalizationFactor void SetDebug bool aDebugValue bool GetDebug When a MatchFilter object is created the template is given and the normalization is done taking into account the two sided PSD Note that the match filtering is performed using zero padding technics The size factor which is given corresponds to the size of the input data which will be filtered The method SetPSD
76. double aDeltaFrequency bool aAlreadyInitialized int LowPassToHighPassTransform The constructor specifies the set of poles and zeros that is the couple Q Q following the convention defined in section 2 5 1 One has also to give the global gain and the sampling frequency of the discrete time sequence to be filtered The SetGain method allows to change the value of the global gain of the filter It also possible to apply the change of gain with a ramp time given in sampling number The GetCurrentGain method returns the value of the global gain The SetGainAtFrequency method allows to fix the gain of the filter at a given frequency The GetGainAtFrequency method returns the gain at a given frequency The ResetHistory method allows to reset all the memory of the filter The SetHistory method allows to set the history of the filter that is the 3 last values of input and the 2 last values of filter output not yet implemented 42 The GetHistory method returns the history of the filter not yet implemented The PrintCoefficients method print out the expression of all the second order filters defined and used in cascade for implemented the requested filter as follow multi pole and zero filter described in section 2 5 2 hhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhh hi Filter composed by 3 filters of second order Global Gain 3 gt Filter No 0 Filter in Laplace space 1 5 30516e 09 s
77. dow PWV computes the pseudo Wigner Ville transform of the windowed data It is the windowed version of the Wigner Ville transform 00 i WV t v e WT a t 7 ea t 7 dr PWV t v gt e 27k a t k xa tk k 3 The only way to obtain a high time resolution is to have aTransformTStep 1 thus incurring a high computational cost It does not seem a recursive algorithm would be possible since the product between data elements are all different from a transform to another For a given transform it is this folding scheme that provides the time resolution whereas frequency resolution relies on window length PWVgauss is personal test of a probably useless modified version of PWV where the window has a Gaus sian profile so as to reduce the terms away from the window centre 138 ST regroups all work done on the S transform and variations internal Welch windowing and non linear template spacing ATrous is a shift invariant version 29 of Dyadic which is the standard discrete wavelet transform DWT though computed using a lifting scheme 30 Both use a dyadic scale meaning the distance between correlated samples is always a power of 2 Starting at the sampling frequency the frequency range is therefore not regularly estimated Better resolution would be possible using the continuous version of the wavelet transform adapted to discrete data The maximum scale is Smaz logs datalength 1 making it necessary
78. e Two AddPix methods are provided to allow for easy addition of pixels Coordinates and values are needed or a Pix object Merging of clusters requires the absorbed cluster to give a pointer to the absorbing cluster using GetAn chor The absorption is called by Absorb that adds all elements from the Piz list of the absorbed cluster to the others before emptying the first list Element extraction is possible with GetPix to retrieve a Pix according to its entry order in the cluster for lack of better discrimination and GetGeoBary that returns the geometrical barycentre of the clusters members Further parameter extraction is possible and straightforward for the most part The GetDim1Spread and GetDim2Spread methods return the maximum coordinate difference between cluster members The GetM11 GetM02 and GetM20 methods provide geometrical moments of second order for the coordinates N 1 mi D Far lt z gt yk lt y gt k 1 140 As for GetEqRadius tr lt gt z M IM mao y lt y gt z Mo2 gt ll m it estimates the radius in number of pixels of a cluster having a Gaussian energy profile highest valued elements at the centre with _13x0 LES 1 so that a one pixel cluster returns 1 Since clusters are built from map pixels a transposition to physical coordinates time and frequency must be done at some point That is the purpose of the SetPhysical meth
79. e Mere pclae Us ena EE A Sy De eeu RE 144 Buste a e at Ry Sy Bye BR Ragin ie Boke TA EE NUE 146 Bilisets tok tedeted tee elle A Pecan tire do esas te eats 146 ASP CTV Ke Les es ae Sa Se Be ee eee he ae pe ee Be Re Se eB 147 ADPOndix 2 A id a nb ashen ER BON Ree Soe SG 149 Appendix fet AR ke a ee A wl el en A a QUE ee WS min 152 Chapter 1 Introduction 1 1 Generalities This library contains several packages dedicated to the search of burst gravitational waves GW It contains some facilities to generate the GW wave forms and place the templates and the algorithms dedicated to look for GW In addition to this burst specific library BuL provides some applications related to the study of GW sources such as spectral estimation noise model generators linear time domain filtering and statiscal algorithms This library is written in C and uses the C standard library and STL classes complex vector string BuL uses also some external packages for example the FFTW C library package version 3 0 1 It is developped on DEC OSF1 V5 2 Linux RH 7 2 It has been successfully compiled with gcc version 3 3 2 All the packages are managed and built using CMT The different packages are the following BuM contains some mathematical functions and some usefull C classes a C class to embed the FFTW C functions for instance BuS contains some C classes to generate some determinist wave forms for Burst events and th
80. e Zwerger and Mueller waveforms catalogue BuF contains a C virtual class to code the different burst search algorithms It contains also an example implementation of the BuD C virtual class and all the different algorithms dedicated to Burst sources ALF PC MF NF BuT contains some C classes to place templates for Bursts search 1D and 2D templates place ment BuC contains some C classes to perform Burst ITF coincidence search Algorithm contains some C virtual classes dedicated to all kind of algorithms except the time doamin linear filters contained in the Filter package LinearFilter contains C classes dedicated to discrte time domain filters Spectral contains C classes dedicated to spectra estimation methods Noise contains some C classes to generate noise white and colored noise ITF contains some header files which include all the characteristics of the dectectors TeF contains C classes dedicated to time frequency algorithms developed for the burst search e Tiler contains C classes providing a geometrical method to tile a 2D parameter space It is used for Damped Sine and Gaussina Sine signals and has been tested for Coalescing Binaries e BuLSet logical package whose role is to manage easily the other BuL packages see Installation section To use the BuL it is mandatory to use the namespace _BuL in all the user s applications 1 2 BuL library installation Since all the
81. e _BuL int main int arg char argc int i j int sampling 20000 int size sampling 2 10 int mean 10 Vector lt double gt data mean size Vector lt double gt out size 2 79 string ar_coeff ar_coeff data AR_Virgo_1024_4Hz dat GaussianAR gene ar_coeff gene SetDebug false gene GetData data gene GetData data for i 0 i lt 10 i cout lt lt data i lt lt endl data WriteFile data dat Spectral spec size Taper Haming mean Spectral QneSided spec ASDWOSA data out sampling spec WriteSpectrumFile out psd dat sampling return 0Q 80 2 7 3 2 Class GaussianARMA Include file Noise_GaussianARMA_hxx Constructors GaussianARMA string aFileName GaussianARMA ITF Detector aITF Public methods int GetData Vector lt double gt amp aData int GetPSDFromData double aDeltaFrequency int aAverageNumber Vector lt double gt amp aPSD int GetPSDFromCoefficients double aDeltaFrequency Vector lt double gt amp aPSD int GetSamplingFrequency int GetPolesNumber double GetSigma void SetSeed unsigned long aSeed void SetSeedLocalTime unsigned long aSeed void SetDebug bool aValue bool GetDebug This class provides a facility to generate a Gaussian Noise whose PSD is parametrized by a ARMA filter model 13 As for the GaussianAR class two constructors are provided for the first one I the user must give the file of the ARMA
82. e 56 2 5 8 Glass GoertZel cs n este a wes de oh BO ad RE TE ES 59 2 6 ITF detector characteristics 62 26 1 Class TLE ini ear ar este te ad BR A he teed Ba ae ee ee ees GE 63 26 27 Class Vit gee enw al a E doe eee Se Doe a Sarah 65 263 Class Hanford mananan Sy EDR MR Qi eh Boe Haw ee ak ded DAs A 67 2 64 Glass Livingston 64 kok Skt ae aida hed alk Re A Aantal ao ee BA ae ture 69 26 0 Clas GEO yi acs Se Se de ea Nee ee ke ae wae re RE Pea hc Ke Se Be ete He ED 70 26 6 Class LAMA 05 66 6 68 6 Ges A eG ee RE RS 71 2680 Class AIGO LA at ieee Fy Paiute ies ane be be Gea bre Se Aer ar aoe a 72 2 7 Noise generators of white and colored noise 73 2 01 Random generators sco ed Gee ee ee ee 73 2 7 1 1 Class Random NR 73 2 7 1 2 Class RandomMT 74 2 7 2 Gaussian and white generators 75 2 7 2 1 Class Gaussian NR 75 2 202 Class Gaussian ME sete 5 8 a RA Doe Die Ae He Ss M EY 76 2 7 3 Gaussian and colored generators 78 2 7 3 1 Class GaussianAR e 78 2 7 3 2 Class GaussianARMA 81 2 7 3 3 Class GaussianColored 82 ZTA COTAS AR is E A hes ap ne en tee re A AAA Ad Hs EE 83 251 07 Class ARMA soe 22d ue d ai A A A pa 85 2 16 Class Gauss
83. e FCUT_QUARTZ_YAG 500 VIRGO geometry angles are in degres static const double LONGITUDE_VIRGO 10 5 static const double LATITUDE_VIRGO 43 63 static const double AZIMUTHAL_VIRGO 206 5 static const double ARM_ANGLE_VIRGO 90 0 Hanford 16384 kHz sampling static const int FZERO_HANFORD 16384 Hanford noise SIGMA_NOISE SENS_NOISE sqrt FZERO 2 static const double SIGMA_NOISE_HANFORD 4 e 21 static const double SENS_NOISE_HANFORD 4 e 23 Hanford_2k expected sensitivity curv parameters static const double 1_2K 1 94e40 static const double S2_2K 1 2e 36 static const double S3_2K 4 68e 46 static const double S4_2K 2 88e 46 static const double FCUT_2K 150 Hanford_4k expected sensitivity curv parameters static const double S1_4K 1 94e40 149 static const double S2_4K 1 2e 36 static const double S3_4K 4 68e 46 static const double S4_4K 2 88e 46 static const double FCUT_4K 150 Hanford geometry angles are in degres static const double LONGITUDE_HANFORD 119 41 static const double LATITUDE_HANFORD 46 45 static const double AZIMUTHAL_HANFORD 261 8 static const double ARM_ANGLE_HANFORD 90 0 Livingston 16384 kHz sampling static const int FZERO_LIVINGSTON 16384 LIVINGSTON noise SIGMA_NOISE SENS_NOISE sqrt FZERO 2 static const double SIGMA_NOISE_LIVINGSTON 4 e 21 static const double SENS_NOISE_LIVINGSTON 4 e
84. e power spectral density which is foreseen or measured for each interferometers The name of the detector is provided by an ITF_Detector enumeration type variable declared in the ITF class and whose definition is given in 2 15 and below Detector Virgo Virgo_Quartz Virgo_Quartz_Yag Hanford_2k Hanford_4k Livingston GEO TAMA AIGO The power spectrum density function one sided PSD of Virgo and TAMA interferometer is parametrized as follow E p S pendulum dh mirror S atace 1 f 5 F hot res The power spectrum density function of LIGO same for the three interferometres and GEO is parametrized as in 16 Note that the AIGO class does not provide any PSD 62 2 6 1 Class ITF Include file FITEFITF hxx Constructors ITF ITF ITF Detector aName Public methods void SetSamplingFrequency int aValue int GetSamplingFrequency int GetOneSidedASD Vector lt double gt aData int GetTwoSidedASD Vector lt double gt aData int GetOneSidedPSD Vector lt double gt aData int GetTwoSidedPSD Vector lt double gt aData int GetOneSidedASD Vector lt double gt amp aData double aFactor int Get woSidedASD Vector lt double gt amp aData double aFactor int GetOneSidedPSD Vector lt double gt amp aData double aFactor int GetTwoSidedPSD Vector lt double gt amp aData double aFactor void PrintOutGeometry double SigmaOfOneSidedPSD Vector lt double gt aPSD double SigmaOfTwoSidedPSD Vector lt double gt aPSD
85. e segment 0 1 Note that contrary to the previous functions there is a change of a new sequence seed at each call of the function For that purpose it uses the SetSeedLocalTime method The GetGaussianData method returns a Gaussian random variable whose value lays in the segment 0 1 The Gaussian variable is built from the random variable using the Box Muller algorithm described in 11 Note that there is no change of a new sequence seed at each call of the function 74 2 7 2 Gaussian and white generators The two following classes allow to generate a white and gaussian distribution noise By default the mean of the distribution is zero and the sigma n is derived from the Virgo performance values stored in the header file ITF_Constant hxx located in the BuL ITF package one assumes that the spectral power density of Virgo is equal to its minimal value h fs 20 kHz sampling frequency 2 Sr f min hn Sa minimal value of the one sided power spectral density hn 4 107 VHz VIRGO minimal sensitivity noise On hay 4 e VIRGO noise sigma for bursts 2 7 2 1 Class Gaussian_NR Include file Noise_Gaussian_NR Constructors Gaussian_NR I Gaussian NR double aSigma II Public methods in addition to the one of the class Random_NR int GetData Vector lt double gt amp aData int GetData Vector lt float gt amp aData double Get Value int SetSamplingFrequency int aValue int GetSamplingFrequ
86. egmentSize x aSegmentNumber Finally the power spectrum is normalized such that the result is in unit Hz and the amplitude spectrum is in unit V Hz The Periodogram method computes simply the square modulus of the FT of the data vector no data windowing no average nor normalization are performed The PSD and ASD methods average the power spectrum over distinct data segments and make use of tapered data while the PSDWOSA and ASDWOSA methods average the power spectrum over 50 overlapping data segments The SpectrumPhase method computes the phase of the data vector When the option average has been chosen in the constructor the SpectrumPhase method averages the FFT and then computes the phase The WriteSpectrumFile methods write out in a file the spectrum given in the arguments associated to the corresponding frequency One can also write out part of the spectrum corresponding to a given frequency band The FlatnessOfPSD method computes the parameter which estimates the spectral flatness of a PSD given as argument 29 Aa ln PSD N af 72 PSD f df If PSD f is peaky 0 If PSD f is flat 1 The SigmaOfPSD method computes the o of the density spectral as following fo 2 o if PSD f df for a two sided PSD fo 2 o 2 2 Gs f PSD f df for a one sided PSD 0 The GetEstimated method returns true if a class to PSD PSDWOSA ASD or AS
87. ency int SetGaussianSigma double aSigma The class Gaussian_NR is making use of the Random_NR class to generate random numbers The first constructor initializes the sigma of the white and Gaussian noise to the official Virgo value The second constructor allows to specify the sigma value The GetData methods fill up a vector of the given size of values generated according to a white noise random distribution Note that the memory allocation of the vector has to be done by the user This method makes use of the Random_NR GetGaussianData Before calling this method GetData changes the seed using SetSeedLocalTime in order to avoid problem The GetValue method returns a scalar value generated according to a white noise random distribution The condition are identical to GetData method The SetSamplingFrequency method changes the sampling frequency private variable and also the sigma value of the noise according to the formula given above The sampling frequency is used only in the computation of the sigma when one uses the constructor I The GetSamplingFrequency method gets the value of the sampling frequency used On the contrary the SetGaussianSigma method changes only the value of the sigma of the noise 75 2 7 2 2 Class Gaussian MT Include file Noise_Gaussian MT Constructors Gaussian_MT I Gaussian MT double aSigma II Public methods in addition to the one of the class Random MT int Ge
88. equired objects and allo cates memory It is mandatory to call the parameterizing method before using any method Repetitive calls of a transform on chunks of data of the same size are possible Analytic computes the analytic signal related to the input data It is the inverse FFT of the data FFT after setting all negative frequencies coefficients to zero It is used in the PWV transform to reduce border effects in frequency When a transform relies on a specific time window shape it is defined by calling either Gaussian Window for Gaussian profiles of sigma TO in seconds centered on the middle of the window or Hanning Window of the shape k a 1 a cos 2 xm 7 for other purposes The window coefficients are not affected when calling the transform therefore redef inition is required only to change their parameters The standard parameters aTransformTLength and aTransformTStep set the size in samples of the win dow applied to the data and the number of samples between consecutive computations of the transform Gabor computes the square of the FFT of windowed data of the set length with the previously defined window profile 00 gabor t v e it t 7 dr O0 A trade off can be achieved between aTransformT Step and the redundancy of the transform and the narrowness of the Gaussian profile The limit for a flat profile is a spectrogram for which the delay between two computations should be the length of the win
89. er gt amp aClusterMainList vector lt Cluster gt amp aClusterList PowerTrigger vector lt Cluster gt amp aClusterList const double aThreshold const string aF PowerAndPeakTrigger list lt Cluster gt amp aClusterMainList vector lt Cluster gt amp aClusterListIn PowerAndPeakTriggerNLSpacing list lt Cluster gt amp aClusterMainList vector lt Cluster gt amp aClus GetAllParameters vector lt Cluster gt amp aClusterList const string aFileName GetAllParameters list lt Cluster gt amp aClusterList const string aFileName GetSomeParameters list lt Cluster gt amp aClusterList const string aFileName GetStat const string aFileName void SetDebug const bool aValue bool GetDebug The DimXSet methods inform of the physical values associated to the map indexes They are used to call the Cluster ToPhysical methods before applying the thresholds PowerTrigger only considers the total energy of a cluster while PowerAndPeakTrigger also checks the maximum pixel energy of a cluster GetSomeParameters and it GetAllParameters only defers by the list of parameters returned but none gives the full available list 145 2 13 BuC to be written 2 14 BuLSet This package is a logical package whose role is to manage and controle the installation of the BuL library It also contains this document 146 2 15 Appendix 1 One has defined the following enum e BuS_Source hxx typedef enum WhiteNoise C
90. er has at least a pole 40 which is not canceled by a zero The frequency response of the filter is given by the transfer function evaluated in the unit cercle This is also the DFT of the impulse response Considering z e the frequency response is k M a H e2 bo k 1 1 Che wA ao TI 1 deriva where A A Ck and dx are complex variables The most general form o fthe transfer function of an IIR with real coefficients is k M k M H z A k 1 1 fez 1 k 1 1 grz D A g z 1 IEM 1 272 REA dez 1 1 diz where M M M and N N Na The first order factors represent real zeros at f and real poles at ck The seond order factors represent the complex conjugate pairs of zeros and poles This expression suggests to implement an IIR filter as a cascade of first order and second order filters One of the advantage is to provide a very general implementation to any order filter provided the second order filter implementation The second advantage is to avoid numerical problem for high order filter design when using the bilinear transform This filter structure has been chosen for the IIR class implementation The main issue in designing a digital ITR or FIR filter is to translate the continuous time analog filter design into the discrete time digital domain That is to say to determine the a and b coefficients from poles and zeros given in the complex plane Let s start with t
91. er map edges or avoid map over writing ClusterCoincidence tests the simultaneity of clusters using strict conditions the cluster with higher central frequency must have start and end times included within the time limits of the lower central frequency one Tentative methods for coherent analysis are included StackClusterEnergy adds to a given list of clusters the energy found at the corresponding pixels of the input map CoherenceHK stacks maps before applying the HK clustering algorithm There is a reference map and possibility of one or many delayed maps For these two consecutive maps are provided in the hypothesis that while map size is the same for all channels delays imply different time indexes will be requested GetClusterList returns the vector of extracted clusters Note that this list can concern the map just processed or the previous one if a connected clustering method was called In the case the latest map s clusters are accessed by GetLastClusterList To erase the clusters lists use FlushMapClusters To punctually suppress cluster map bridging call FlushMap Connections 2 12 5 Class Statistics This class is in charge of both thresholding of cluster parameters and extraction of the values for the retained cluster and their copy to file At the moment the methods extracts a fixed set of parameters The parameters in use at the moment are the following e Cluster size e Cluster weight the sum of pixel values in
92. er to get a normal distribution o 1 in the hypothesis of a gaussian noise Practically the templates are normalized PCi t The Threshold method compares the filter output vector to the threshold value and then produce a list of micro events if to be This list of micro events can be clusterized if necessary to reduce the false alarm rate It is recommended to use one of the algorithm devoted to Correlator filters ClusterCorrela torLoose or ClusterCorrelatorTight methods to perform the event clusterization The method EnableRefreshPSD method allows to disable or enable the estimation of the PSD which is done each time one calls the Amplitude or detect methods The method UpdateCurrentPSD method allows to update the internally defined PSD at any moment giving a vector of data for the estimation The method SetThreshold changes the initial values of the thresholds defined for each window The method GetThreshold returns the value of the thresholds The method GetNWindow returns the number of windows which are used in paralell The method GetWindow int alndex returns the size of the window corresponding to the index alndez 2 11 9 Class DSC to be written 2 11 10 Class PF to be written 2 11 11 Class GDF not yet available 135 2 12 TeF 2 12 1 Intro Time frequency transforms all propose to layout local signal frequencies We can already distinguish between sh
93. et the GW amplitude corresponding to the 78 ZM simulated waveforms using equation 2 25 and assuming we are located in th equatorail plane 0 5 The waveforms of the ZM catalogue 18 have been downsampled at 20 kHz The list of the 78 waveform file names is given in 2 17 for all waveforms excecpt the one with 8 0 4 which are not considered in the 78 ZM waveforms 97 Include file BuS_ZM hxx Constructors ZM string aName double aSNR ZM string aName int aDistance ZM string aName double aSNR int aSampling Vector lt double gt amp aPSD ZM string aName int aDistance int aSampling Vector lt double gt amp aPSD ZM string aName double aSNR ITF Detector alTFName ZM string aName int aDistance ITF Detector alTFName Public methods int GetData double aPosition Vector lt double gt amp aData int GetData double aPosition Vector lt float gt amp aData int SetZM string aName int SetDistance int aDistance The method GetData fills the data vector with the chosen ZM signal which is located at the relative position given by the argument aPosition This position is a fraction of the data vector size If the size of the ZM signal is longer than the vector data size a warning message is printed out and the vector is not filled up The method SetZM changes the waveform by a new one The method SetDistance initializes the distance and the normalization factor is recomputed Program exam
94. f the 25 DFM waveforms 152 Bibliography 10 11 12 13 14 15 16 17 18 19 20 CMT v1r9 C Arnault Configuration Management Tool http www lal in2p3 fr SI CMT htm Doxygen M Frigo S G Johnson FFTW v2r13 The Fastest Fourier Transform in the West User s Manual available at http www fftw org http tycho usno navy mil time html http hpiers obspm fr eop pe D B Percival A T Walden Spectral analysis for physical applications Cambridge University Press 1993 E Cuoco A Vicer Mathematical Tools for Noise Analysis Virgo DAD 1999 P D Welch IEEE Transactions on Audio and Electroacoustics 15 70 3 1967 E Chassande Mottin Testing the normality of the gravitational wave data with a low cost recursive estimate of the Kurtosis gr qc 0212010 A V Oppenheim R W Schafer Discrete time signal processing Prentice Hall Signal Processing Series New Jersay 1999 L R Rabiner B Gold Theory and application of digital signal processing Prentice Hall New Jersay 1975 W H Press S A Teukolsky W T Vetterling B P Flannery Numerical Recipies in Fortran Cambridge University Press 1992 http www math keio ac jp matumoto emt html C W Therrien Discrete Random Signals and Statistical Signal Processing Prentice Hall New Jersey 1992 E Cuoco G Curci Modelling a VIRGO like noise spectrum VIR NOT PIS 1390
95. f the characteristic relaxing time of the template signal The Threshold method provides at the end a empty if not found list of micro events for all the windows 122 It has been observed that when an true Burst event is present in the data it often produces a set of micro events for several which tend to cluster In other word these micro events seem to correspond to the same physical event although this behavior is barely attented in the case of a false alarm only one window trigger and of very short duration for example If one clusters the micro events one can reduce easily the number of false alarm by comparing the size of the events which are supposed to be larger in the case of a physical event compared to a false alarm Note that the clusterization of the micro events will depend on the filters family basic correlator and will be performed by the user using the methods proposed in the class Event question en suspens doit on faire une classe Cluster ou bien les algos de clusterization sont ils des methods non membres de la classe Event ou Detector The Detect method provides a direct way to find out the events the micro events clusterization is also performed and produce a list of events giving the vector of data as input argument Note that the Detect method returns also the array of amplitude Finaly the End method is foreseen to provide some final infor
96. g the corresponding methods see TilerWithEllipses section 109 2 10 2 Class TilerPoint A simple class to define points in a n dimensional space It is used by other Tiler classes to define center of patterns or to test the parameter space Include file TilerPoint hxx Constructors TilerPoint int aDim 2 TilerPoint const TilerPoint amp aPoint Destructors TilerPoint Public methods double GetCoord int aN int SetCoord int aN double aValue int GetDim static int GetMaxDim double Distance TilerPoint amp aPoint The constructor TilerPoint int aDim 2 creates a point in a aDim dimension space by default set to 2 and set it at the origin The constructor TilerPoint const TilerPoint amp aPoint creates a point which is an exact copy of aPoint dimension and coordinates The GetCoord int aN method returns the aN coordinate It returns TilerPoint PointError if aN is out of range The SetCoord int aN double aValue method sets the aN coordinate to a Value The GetDim method returns the dimension of the space The GetMaxDim method returns the maximal allowed dimension for the space The Distance TilerPoint amp aPoint method computes the euclidian distance between the current point and aPoint 110 2 10 3 Class TilerPattern This is a generic class internally used by Tiler for the definition of geometrical patterns No use is foreseen for standard users of the T
97. gt amp aPSD int aSampling int aFMin int aFMax double Sigma0fPSD Vector lt double gt amp aPSD int aSampling bool GetEstimated void SetDebug bool aValue bool GetDebug Several options are proposed to estimate a PSD e one sided or two sided PSD e with or without taper e with or without averaging 28 e with or without Welch averaging approach WOSA The one sided and the two sided PSD are normalized such as N aSegmentSize e Two Sided _ lt FT f gt _ fo N fo fo PSD f m He f Big DR UN e One sided e EAS fo fo PSD f 2x ET fa Beng PSD 0 oti f 0 0 considering that lt FT f gt is the average over one or several FT data segments Constructor arguments e aSegmentSize size of the data segment which undergoes the Fourier Transform The Two sided PSD output vector size is equal to aSegmentSize although the one sided PSD output vector size is aSegmentSize 2 e aWindow taper window type see class Taper for the list of available windows e aSegmentNumber number of sub segments of the input data vector which will undergo a FT and then be averaged In the case of the PSD estimation with the Welch method the effective number of sub segments is 2 x aSegmentNumber 1 e aPSDType one sided Spectral OneSided or two sided Spectral TwoSided PSD Important remarq Be aware that the size of the data vector which is given as input of the different methods has to be equal to the product aS
98. gt amp aWindow int GetRawDataWelch Vector lt double gt amp aWindow int GetRawDataBartlett Vector lt double gt amp aWindow void SetDebug bool aDebug Value bool GetDebug The class Taper defines a enumeration type variable Window whose definition is given below and in section 2 15 Window None Hann Welch Haming Bartlett The constructor Taper creates the object without defining the size of the window This means that the size must be assigned later during a object copy for example The other constructor Taper aSize creates the Taper object defining the size of the window The methods GetDataHann fill the vector with the appropriate taper defintion The window is normalized such as having w 1 The methods GetRawDataHann fill the vector with the appropriate taper defintion The window is not normalized Program exemple include math h include Spectral_Taper hxx using namespace _BuL int main int narg char argc int i int size 128 Vector lt double gt a size Vector lt double gt b size Taper taping size 26 taping GetDataHann a taping GetDataHaming b a writeFile hann dat b writeFile haming dat return 0 27 2 3 3 Class Spectral This class provides several methods to perform spectral estimation All the methods are based on the Fourier Transform and hence make a use of the FFT class The main interest of the Spectral object
99. h gt include lt math h gt include lt fstream gt include DCPP_Vector hxx include LinearFilter_Butterworth hxx include Spectral_Taper hxx include Spectral_Spectral hxx 53 include Noise_GaussianAR cxx using namespace _BuL Test program for pass band filtering int main int argc char argv double tf_pass tf_stop_low tf_stop_high frequency_pass_low frequency_stop_low double frequency_pass_high frequency_stop_high double sampling 20000 int datasize 20000 Vector lt double gt module Vector lt double gt phase Vector lt double gt data datasize Vector lt double gt output datasize Vector lt double gt out datasize 2 tf_pass 99 tf_stop_low 01 tf_stop_high 01 frequency_stop_low 100 frequency_pass_low 200 frequency_pass_high 300 frequency_stop_high 405 Spectral psd datasize Taper Hann 1 Spectral 0neSided Butterworth filter sampling frequency_stop_low frequency_pass_low frequency_pass_high frequency_stop_high tf_stop_low tf_pass tf_stop_high filter PrintCoefficients filter GetTransferFunction module phase 1 module WriteFile mod dat phase WriteFile phase dat Gaussian colored noise char Envvar Envvar getenv NOISEROOT string Ar_coef Envvar Ar_coef Ar_coef data AR_Virgo_1024_4Hz dat 54 GaussianAR gene Ar_coef gene GetData data data WriteFile data
100. he MA coeffients The GetPSDFromData method estimates the onde sided PSD of a vector of generated data using the ARMA parameters of the class The GetPSDFromCoefficients method computes the one sided PSD numerically using the ARMA parameters following ae bye 2Mifk P L Di ane np The GetSamplingFrequency returns the Sampling frequency used in the AR process PSD f 2 x Sigma x The GetPolesNumber returns the number of poles number of AR coefficients of the AR process The GetSigma returns the value of the discrete sigma which is used in the AR process The SetSeed allows to set the internal seed number to the given value This method can be called at any instant The SetSeedLocalTime allows to set the seed number used by MT algorithm of a given number which depends on the local time with a precision of 1 s The internal seed number is set to the value given as argument This method can be called at any instant 2 7 3 3 Class GaussianColored Include file Noise GaussianColored hxx Constructors GaussianColored int aSize ITF Detector alTF Public methods int GetData Vector lt double gt amp aData int GetData Vector lt float gt amp aData int SetSamplingFrequency int aValue int GetSamplingFrequency void WriteFileGeneratedASD string aName void WriteFileGeneratedASD int aSize string aName double GetValue 82 2 7 4 Class AR Include file Noise_AR hxx Constructors A
101. he value of the sampling frequency The value of the sigma of the white noise is also changed according to the formula given in section 2 7 2 The normalization factor in the case of SNR normalization is recomputed The method SetSigmaNoise method changes the value of the white noise sigma in the case of a white noise option choice The normalization factor in the case of SNR normalization is recomputed The method GetSigmaNoise method returns the value of the sigma of the white noise The following methods are active only in the case of the colored noise option The method SetPSD changes the definition of the PSD and the sampling frequency The method SetITF changes the definition of the interferometer whose noise characteristics are taken into account for the source normalization The following methods are active in the case of the both noise options The method GetSamplingFrequency returns the sampling frequency The method SetSNR changes the value of the SNR The normalization factor is then recom puted 90 2 8 2 Class Gaussian The class Gaussian provides facilities to generate a waveform 12 f t ae Include file BuS_Gaussian hxx Constructors Gaussian double aSigma double aSNR Gaussian double aSigma double aSNR int aSampling Vector lt double gt amp aPSD Gaussian double aSigma double aSNR ITF Detector aName Public methods int GetData double aPosition Vect
102. he z domain transfer function ao a127 a227 A z 2 2 1 biz baz 5 There exists several methods to derive the coefficients In the IIR class we have used the bilinear transform which maps region of s complex space to region of z complex plane the left half s plane is mapping the interior of the z unit circle The pure imaginary s axis is mapped onto the z unit circle sA 1 2 2 14271 _ 1 As 2 1 As 2 A Bs Cs H s t D Es Fs 36 The biquadric filter coefficients are C A B1 E ET A C 1 mel C A B1 Sl tg gale D Ful n ly xl E 1 mele ale 1 F E Besides one defines second order filter in s domain as follow 2 8 8 late 2 6 The roots are s 3911 F y 1 4Q The main interest of this notation is that it ensures the stability of the filter since the real part of the poles are always negative We have adopted the following convention to describe a second order and a first order filter using equation 2 6 e second order 040 and Q 0 2 1 22 A 1 _ 1 50 1 Ges Q2 e first order Q 40 and Q 0 1 pole zero at Q A 1 1 B Q C 0 e first order Q 0 s pole zero at 0 A 0 B 1 C 0 37 e first order Q 1 1 pole zero at oo Ques Il O m The last point which is worth mentioning concerning the bilinear transform technics is the frequency warping issue 9 The relation between the continuous time frequency 2 s o iQ and the discrete
103. ian Line ci a e ee 87 201 Class ie id iaa dde dopo lo de Lo Li Ade oa 87 2 7 6 2 Class GaussianLine s osos e e a 88 2 8 BuS sources of Burst events 89 281 Class SQUICE sos Sar Seo sei oe nd a ae rap ae bel gad ged ee ee WE US a E We o ah 90 A Ol CEE CE 91 2 8 3 Class DampedSine 93 2 8 4 Class GaussianSine 95 21820 Class AM nes Mis in Bich A Dee ea de acc ei ve at ale how la 97 2965 Glass DENT ut ii doth GG BE BE Dpt ee SU aa ee 4 99 2 8 7 Class Sim re E epee ede ae tek EE OR OR we Ee GE ae ee e 101 2 9 BuT burst template generators 102 2 9 1 Class E Gaussian 2 205 jae vive Goatees ANS ep Se WE Red aoe eee ee 102 2 92 Class TDamp dSine s 6 0 5 608 unes a ey OY Wow BT ee 106 2 10 Tiler Placement of templates in the parameter space 107 210A Introduction LAA A A a TA A ANN cid Re 107 2101 1 Tiling for malito Aa Ae 107 2 10 1 2 The tiling problem 108 2 10 1 3 Tiler Algorithm 108 210 2 Class TilerPoint viii A Be Ree OE dei an a En Re os 110 2 10 3 Class Lil rPatt tn sai a ge EE oe noue LR ees 111 2 10 4 Class TilerEllipse 4 112 2 10 5 Class TilerRectangle 114 210 6 Glass Pier ns ce fe ee nn AMD
104. ibed in the section of the ITF class Program exemple include ITF_TAMA hxx using namespace _BuL int main int arg char argc int i j int size 10000 Vector lt double gt psd_sdt size TAMA std 10000 std GetOneSidedASD size psd_sdt for i 0 i lt size 1000 i cout lt lt vector ASD lt lt i lt lt lt lt psd_sdt i lt lt endl return Q 71 2 6 7 Class AIGO Include file TTF_AIGO hxx Constructors AIGO AIGO int aSamplingFrequency Public methods void SetSamplingFrequency int aValue int GetSamplingFrequency double GetLatitude double GetLongitude double GetAzimuthal double GetArmAngle void SetLatitude double aValue void SetLongitude double aValue void SetAzimuthal double aValue void SetArmAngle double aValue void SetDebug bool aValue bool GetDebug The AIGO constructor consider the ITF AIGO detector with its nominal fo xxx Hz sampling fre quency The other constructor allow to choose a different sampling frequency if needed The methods are described in the section of the ITF class Program exemple include lt string gt include ITF_AIGO hxx using namespace _BuL int main int arg char argc int i j AIGO std 10000 double longitude std GetLongitude double latitude std GetLatitude double azimuthal std GetAzimuthal double arm_angle std GetArmAngle cout lt lt AIGO Longitude lt lt longitude
105. iler package Include file TilerPattern hxx Constructors TilerPattern int aDim 2 Destructors TilerPattern Public methods void SetCenter TilerPoint amp aPoint TilerPoint amp GetCenter bool Characteristic TilerPoint aPoint void SetCharact bool aFunction TilerPoint amp aPoint double Area The constructor TilerPattern int aDim 2 creates a default pattern in a aDim dimension space by default set to 2 The default pattern is an hypercube of side equal to 1 centered on the origin The SetCenter TilerPoint amp aPoint method allows to define the center of the pattern The center is set to aPoint The GetCenter method returns the center of the pattern The Characteristic TilerPoint amp aPoint method allows to test if aPoint belongs or not to the pattern This function returns true is aPoint belongs to the pattern and false if not The SetCharact bool aFunction TilerPoint amp aPoint method allows to define the characteristic function of the pattern The Area method returns the area of the pattern 111 2 10 4 Class TilerEllipse This class inherits from TilerPattern It defines the elliptic pattern As TilerPattern it is only for internal use of the Tiler package No use should be done by standard users Include file TilerEllipse hxx Constructors TilerEllipse double aX2Coeff 1 double aY2Coeff 1 double aXYCoeef 0 TilerEllipse const TilerEll
106. ime since the previous print By default the DecreaseSurfaceFactor is set to 10 The SetDebugLevel int aValue method sets the debug level Higher is the level the more you get information By default the debug level is one The maximal value is 3 The GetDebugLevel method returns the debug level The SetBoundMin int aIndex double aValue method sets the afndex minimal bound of the rectangle containing the parameter space to aValue The GetBoundMin int alndex method returns the alndex minimal bound of the rectangle containing the parameter space The SetBoundMaz int alndex double aValue method sets the afndex maximal bound of the rectangle containing the parameter space to aValue The GetBoundMax int alndex method returns the afndex maximal bound of the rectangle containing the parameter space The SetCharacteristicFunction bool aFunction TilerPoint amp aPoint method sets the charateristic function of the parameter space This function should return true or false if the point aPoint belongs or not to the parameter space The SetMarginCharacteristicFunction bool aFunction TilerPoint amp aPoint method sets the chara teristic function of the margin around the parameter space This function should return true or false if the point aPoint belongs or not to the margin The SetStartingCenter double aCoord method allows to define the coordinate of the starting point of the tiling procedure The aCoord a
107. int WriteOneSidedSpectrumFile Vector lt double gt aPSD string aFileName int Write TwoSidedSpectrumPFile Vector lt double gt aPSD string aFileName double GetLatitude double GetLongitude double Get Azimuthal double GetArmAngle void SetLatitude double aValue void SetLongitude double aValue void SetAzimuthal double aValue void SetArmAngle double aValue void SetDebug bool aValue bool GetDebug The ITF constructor defines an empty object with null geometrical characteristics while the TF ITF Detector aName creates an object filling the sampling information and geometrical characteristics corresponding to interferometer name The SetSamplingFrequency method allows to re define the sampling frequency The GetSamplingFrequency method returns the current sampling frequency The GetOneSidedPSD and GetOneSidedASD methods fill the input data vector of size N with the One sided PSD or ASD The argument factor allows to get a PSD or ASD which is scaled by this factor The frequency resolution is given by Fsampling f ON The GetTwoSidedPSD and GetTwoSidedASD methods fill the input data vector of size N with 63 the Two Sided PSD or ASD The argument factor allows to get a PSD or ASD which is scaled by this factor The frequency resolution is given by f Promoting The SigmaOfOneSidedPSD and SigmaOfTwoSidedPSD methods compute the sigma of the PSD according to 2
108. ions are that the center belongs to the parameter space and it does not fall in an ellipse already 108 Containing Rectangle Figure 2 19 Definitions of the parameter space the margin space and the containing rectangle in the list It is also checked that the ellipse covers a part of P not already covered by other ellipses The minimal surface to be covered is 0 by default and can be changed by the SetMinimalCoveredSurface method This test is made doing a random generation of points in the ellipse the number of tests can be defined by the SetNbLoopsComputeCoveredSurface method By default 10000 tests are performed If the ellipse is useful its center is added to the center list and will be used later to place other centers This iterative algorithm stops either when the full surface of the parameter space is covered or when no new ellipse can be placed any longer A second iteration is performed to treat the edges of the parameter space The same procedure of placement is followed but centers outside the parameter space but inside the margin space see Figure 2 19 are accepted Concerning the placement of the six neighbors three parameters can be adjusted by the user the angle offset o the turn direction TD 1 and the safety factor SF In the unit circle referential the nt neighbor is placed in polar coordinates at SF 9 TD n 7 3 By default 6 0 TD 1 and SF 1 Theses values can be set by the user usin
109. ious 0 Remaining Surface 2 24446 Used Ellipses 11 Total Time Elapsed 0 Since Previous 0 Remaining Surface 1 90807 Used Ellipses 13 118 Total Time Elapsed 0 Since Previous 0 Remaining Surface 1 43244 Used Ellipses 16 Total Time Elapsed 0 Since Previous 0 Remaining Surface 1 13078 Used Ellipses 18 Total Time Elapsed 0 Since Previous 0 Remaining Surface 0 702761 Used Ellipses 21 Total Time Elapsed 0 Since Previous 0 Remaining Surface 0 371225 Used Ellipses 23 Before treating edges surface not covered 0 371225 4 gt 9 28062 Before treating edges Number of ellipses 23 End of Tiling process Surface not covered 0 142635 4 gt 3 56587 Number of ellipses 35 Nb ellipses 35 Space parameter area 4 Area covered by the ellipses in the parameters space 3 85737 Ratio 96 4341 Total area covered by the ellipses 6 87223 Ratio 171 806 Write kumac file Write kumac file Done Centers List To File center txt Write Centers List file Done Figure 2 20 Tiling of a triangle by circles 119 2 11 BuF burst filters This package contains a Detector class which is a virtual class that provides the different methods which allow to define and apply the Burst algorithms to look for events hiden in chunk of data Indeed for all the Burst algorithms it is mandatory first to define the input parameters number of windows or templates values of the threshold then to cal
110. ipse amp aEllipse Destructors TilerEllipse Public methods void SetX2Coeff double aValue void SetY2Coeff double aValue void SetXYCoeff double aValue double GetX2Coeff double GetY2Coeff double GetXYCoeff double GetLongAxis double GetSmallAxis double GetTiltAngle bool Characteristic TilerPoint aPoint double Area An ellipse centered on xo yo is defined in the following way E 1f 2 y ER X2Coef f x zo Y2Coef f y yo 2 x XY Coef f x zo y Yo lt 1 Be careful with the factor 2 in the xy coefficient The constructor TilerEllipse double aX2Coeff 1 double aY2Coeff 1 double aXYCoeef 0 creates an ellipse centered on the origin with X2Coeff aX2Coeff Y2Coeff aY2Coeff and XY coeff aXYCoeef By default an unit circle is created The constructor TilerEllipse const TilerEllipse amp aEllipse creates an ellipse which is an exact copy of aEllipse center and ellipse coefficients The SetX2Coeff double aValue method sets X2Coeff to aValue The SetY2Coeff double aValue method sets Y2Coeff to aValue The SetX YCoeff double aValue method sets XYCoeff to aValue The GetX2Coeff method returns the X2 coefficient The GetY2Coeff method returns the Y2 coefficient The GetX YCoeff method returns the XY coefficient 112 The GetLongAzis method returns the value of the ellipse long axis The GetSmallAzis method retu
111. is estimation algorithm is not optimized since it assumes that the slice step A parameter may be different from 1 which implies no recursive relations between two estimations The Computelterative method compute the Kurtosis quantity using the input Vector and fill the result in the output Vector This Kurtosis estimation algorithm is more optimized than in the previous method since it assumes that the slice step A parameter is equal to 1 which implies some recursive relations between two estimations The ComputeRecursive method compute the Kurtosis quantity using the input Vector and fill the result in the output Vector This fast Kurtosis estimation algorithm implements the method described in 8 This is not yet coded 2 4 2 Mohanty s algorithm 34 2 5 LinearFilter time domain linear filters This package provides different classes for designing and applying digital filter in time domain 2 5 1 Introduction We consider in this package linear digital filters which are linear and time invariant LTI system applied on discrete time sequence z n A LTI system is entirely defined by its impulse response That is the response to an impulse d n k k co k 00 k o00 T ekon Y Ton Y elan k 0oo0 k 0o0 k 0o0 y n is then written as a convolution sum The LTI filters belong to two different categories e FIR Finite Impulse Response the impulse response has a finite number of non zero values This kin
112. ise Program example include Noise_Gaussian_MT hxx include DataContainer_Vector hxx include Spectral_Spectral hxx using namespace _BuL int main int arg char argc int size 2000000 int mean 10 Vector lt double gt data size Vector lt double gt out size 2 mean Gaussian_MT gene new Gaussian_MT gene gt SetGaussianSigma 1 gene gt GetData data 76 data WriteFile white_gaussian dat Spectral spec size mean Taper Haming mean Spectral OneSided spec ASD data out 20000 spec WriteSpectrumFile out psd dat 20000 return Q 77 2 7 3 Gaussian and colored generators The 3 following classes provide different facilities to generate a Gaussian and Colored noise according to different model The Two first classes supposed that the noise can be modelized by a Auto Regressive AR or Auto Regressive Moving Average ARMA models 14 15 13 The last generator performs the convolution of a white and Gaussian noise with a theoritical processus whose PSD is known This method while efficient requests to perform 2 Fourier Transforms time consuming and is not able to produce a continuous chunk of colored and Gaussian data since the phase is not set properly 2 7 3 1 Class GaussianAR Include file Noise_GaussianA R_hxx Constructors GaussianAR string aFileName GaussianAR ITF Detector alTF Public methods int GetData Vector lt double gt amp aData i
113. mac std string aFileName int aPeriod method writes a kumac file named aFileName for PAW The following lines show an example of the kumac contents set faci 37 ellipse 3 71877 0 918781 0 0114446 0 301818 51 906 set faci 38 ellipse 3 731 0 934381 0 0114497 0 301136 51 908 wait set faci 39 ellipse 4 08118 0 64727 0 011079 0 312601 52 6857 set faci 40 ellipse 3 74324 0 949989 0 0114541 0 300423 51 9139 The period of the wait command is given by aPeriod 2 10 7 1 A TilerWithEllipses example The Tiler package contains two examples test_Tiler and Test_Triangle The first one tiles a rectangle of given dimensions with unit circles All circle centers must belong to the rectangle The second test program tiles a squared triangle ABC squared in B with A 2 2 B 2 4 C 6 4 using circles of radius 0 25 No circle center can be put under the AC line Here is an example of the output of the program and Figure 2 20 shows the result in the parameter space csh gt test_Triangle exe 3 3 1 Start Tiling procedure Parameter space Area 4 Total Time Elapsed 0 Since Previous 0 Remaining Surface 3 80365 Used Ellipses 1 Total Time Elapsed 0 Since Previous 0 Remaining Surface 3 44484 Used Ellipses 3 Total Time Elapsed 0 Since Previous 0 Remaining Surface 3 09738 Used Ellipses 5 Total Time Elapsed 0 Since Previous 0 Remaining Surface 2 74301 Used Ellipses 8 Total Time Elapsed 0 Since Prev
114. mation print out or files tuples filling for example about the result of the algorithm applied on a chunk of data The Statistique method writes in a file some information about the list of event 123 2 11 3 Class Event An event is defined by its e starting time Tstart e ending time Tena e maximal amplitude value Amaz e maximal amplitude time Tinas e bin tolerance a e name of the algorithm which has defined the event e type of the algorithm e numero of template which has triggered the event e number of micro events which compose the event e list of micro events which have been merged in the event e GPS time approximation 5 107 s All the quantities which concern time information are sampled in bin and are then relative to the first bin of the chunk of data which is analyzed The bin tolerance field is defined as the window size of the window which has triggered the event for basic filters or as twice the gaussian signal template width for the Peal Correlator filter This information is used by the clusterization algorithm as it is explained in the following Note that when an event is the result of a clusterization of several micro events the bin tolerance is redefined according to the clusterization algorithm as well as the numero of the template and the number of micro events which compose the event is set to the right value Finally we keep an history of the micro events which are
115. mpling itf GetOneSidedPSD PSD on calcule les coeff a partir de la PSD algo ComputeAR PSD PSD_scale_factor f_cut cout lt lt sigma lt lt sqrt algo GetSigma poles_number lt lt endl algo SaveAR AR1 dat on calcule les coeff directement a parti du nom du detecteur algo SetDeltaFrequency 1 algo ComputeAR ITF Virgo f_cut cout lt lt sigma lt lt sqrt algo GetSigma poles_number lt lt endl algo SaveAR AR2 dat To determine the optimla value of the poles_number algo SetDeltaFrequency 1 int opt algo FindOptimal0rderMDL ITF Virgo f_cut 20000 cout lt lt Valeur optimale de 1 ordre selon le MDL critere lt lt opt lt lt endl return 0 84 2 7 5 Class ARMA Exliquer ici la procedure de resolution des equations Normales Include file Noise ARMA hxx Constructors ARMA int aSampling int aPolesNumber int aZerosNumber Public methods int ComputeARMA Vector lt double gt amp aPSD double aPSDScaleFactor double aMinimalFrequency int ComputeARMA ITF Detector aITF double aMinimalFrequency int SaveARMA string aFileName void GetAR Vector lt double gt amp aAR void GetMA Vector lt double gt amp aMA int GetSamplingFrequency int SetSamplingFrequency int aSampling int SetDeltaFrequency double aValue double GetDeltaFrequency int SetPolesNumber int aPolesNumber int SetZerosNumber int aZerosNumber int GetPolesNumber in
116. n a large bank of templates in parallel Then the first problem to be solved is to design this bank in a suitable way The usual prescription is to insure that the minimal match MM minimum taken over the full bank between any signal and the templates is above a value 0 97 is a typical choice 2 10 1 1 Tiling formalism As a linear filtering consists in correlating the detector output s t with the corresponding filtering function y t one defines the following scalar product to represent the filtering operation slo ar f 4 AOR 2 28 ond where Sy is the one sided power spectrum density and the symbol means Fourier transform The normalization is chosen so that in case of signal alone s s is equal to the SNR The ambiguity function between two templates ky et ky a is defined by E X dX hs Iran 2 28 The templates are normalized kx kxy 1 kx ax kx ax T is thus a measurement of the closeness between two templates It can also considered as a way to see how well the template kx can be used to detect the signal kx x the ambiguity function is the mean fraction of the optimal SNR achieved at a given distance in parameter space If dX is small enough the ambiguity function can be approximated by a second order power expansion de ea 1 r X dX 1 59 dd dX 2 30 with guv defining a metric on the parameter space P Finally one defines the minimal match MM as the lower bound of the recovered SNRs
117. n for the five constructors one has to indicate e the size of the input vector e the type of DFT FFT_Real when the input vector is real FFT_Complex when the input vector is complex In case of real input data the FFT_Real choice leads to a gain of factor 2 w r t to FFT_Complex choice When creating a FFT object if the direction FFT_Forward FFT Backward is not specified two plans are created one corresponding to the Forward DFT and the other to the Backward DFT The other optional parameter conncerns the planner wisdom If nothing is specified FFT is using the lower level algorithm FFT_Estimate The other alternative options are FFT_MEASURE FFT_PATIENT and FFT_EXHAUSTIVE Once created a FFT object is used to execute a Forward or a Backward DFT on a given vector of data The type of data can be float or double but the user has to take care of compiling his program with the corresponding version of the FFTW library the binary directories are labelled by lt tag gt float for a float version of FFTW and lt tag gt for the default double version of FFTW library By default all the BuL packages used the double version of FFTW The input of a complex data vector DFT Forward and Backward is stored in a complex data vector according to the following scheme e N is even index 0 FT for f 0 index 1 DFT output for f de Ea N fo index 1 N 1 DFT output for f 4 1 amp
118. nes the histogram either giving the bin number or the constant width of the histogram bins The aWeight argument allows to give a specific weight for each entry usually 1 The aScale argument is used to rescale the value entered using the method Add The method Add aVal adds en weighted entry in the bin corresponding to aVal that means at bin i such as aMinVal x aBinSize lt aVal lt aMinVal i 1 aBinSize The method GetValue aBinIndex returns the value of the bin number aBinIndez The method GetMean returns the mean of the histogram The method GetSigma returns the sigma of the histogram The method GetEntries returns the total number of entries The method GetBinNumber returns the number of bins 18 The method GetUnderflowNumber returns the number of underflow entries The method GetOverflowNumber returns the number of overflow entries The method GetWeight returns the weight factor The method GetScale returns the scale factor The method WriteFile writes in a file the content of the histograms in a row wise format main and overflow and underflow vectors if the switch has been selected The method WriteHistogramFile writes in a file the content of the histogram in a column wise format The method WriteBinHistogramFile writes in a file the bin value and its content of the histogram in a column wise format Program exemple include BuM_Histogram hxx using namesp
119. nt GetSamplingFrequency int GetPolesNumber double GetSigma int GetPSDFromData double aDeltaFrequency int aAverageNumber Vector lt double gt amp aPSD int GetPSDFromCoefficients double aDeltaFrequency Vector lt double gt amp aPSD void SetSeed unsigned long aSeed void SetSeedLocalTime unsigned long aSeed void SetDebug bool aValue bool GetDebug This class provides a facility to generate a Gaussian Noise whose PSD is parametrized by a AR filter model 13 Two constructors are provided for the first one I the user must give the file of the AR coefficients which correspond to the PSD of the data he wants to simulate to know how to determine such coefficients see section 2 7 4 In the second constructor II the user gives the name of the interferometer wich corresponds to the data he wants to generate All the AR coefficients files corresponding to each official interferometer PSD curves have been produced The format of the AR coefficients file is the following A comment line is indicated with Info The Sigma value is the discrete sigma value The Sampling value gives the sampling frequency of the process The Frequency info line gives the minimal frequency which is simulated Below this frequency the spec trum is flat The AR value gives the number P of poles of the AR filter The coefficients values are then given at the end of the file Sigma 1 5900433817564251e 20 Sampling 20000 Info
120. nt aWindow Vector lt double gt amp aOutput The method SetThreshold changes the initial values of the thresholds defined for each window The method GetThreshold returns the value of the thresholds The method GetNWindow returns the number of windows which are used in paralell The method Get Window int alndex returns the size of the window corresponding to the index alndez 130 2 11 6 Class MF The class MF implements the algorithm described in 24 The algorithm is coded in the public method meanCompute This function computes the mean of N points over a moving window of size n The output result is returned in a vector filled up to bin N n 1 The end of the vector from bin N n up to bin N 1 is filled with zeros The filter output follows a gaussian probability distribution The filter output is given by i k n ds MF Note that the output for each window is normalized and the threshold level may be the same for all the windows The test application source file is test BuF_MF cxx Include file BuF_MF hxx Constructors ME string aName Public methods int Configure string aFileName int Configure Vector lt int gt z aWindowSize Vector lt double gt amp aThreshold int aSampling int Amplitude Vector lt double gt z aData Matrix lt double gt z aAmplitude int Threshold Matrix lt double gt amp aAmplitude list lt Event gt amp aEventList int Threshold Matrix lt double
121. oCorrParametrize const int aTransformTLength const int aTransformTStep AutoCorr Vector lt double gt amp aData Matrix lt double gt amp aMap PWVParametrize const int aTransformTLength const int aTransformTStep PWV Vector lt double gt amp aData Matrix lt double gt amp aMap PWVgauss Vector lt double gt amp aData Matrix lt double gt amp aMap STransformSetEdge const int aEdgeSize STransformParametrize const int aMinFreq const int aMaxFreq const double aTransfor STransform Vector lt double gt amp aData Matrix lt double gt amp aMap STransform Vector lt double gt amp aData Matrix lt complex lt double gt gt amp aMap bool STransformProgressive Vector lt double gt amp aData Vector lt double gt amp aSTed 11Qnly dyadic down sampling is available yet 120ne of the reasons why frequency to Hz and time conversion to GPS are not implemented yet 136 int int int int int int int int int STransform Vector lt double gt amp aData Matrix lt double gt amp aMapAbs Matrix lt double gt amp aMapAr STransformWelchParametrize const int aMinFreq const int aMaxFreq const double aTra STransformNLSpacingWelchParametrize const int aMinFreq const int aMaxFreq const do STransformASigmaWelchParametrize const int aMinFreq const int aMaxFreq const doubl STransformWelch Vector lt double gt amp aData Matrix lt double gt amp aMap1 Matrix lt double gt amp aMa ATrousParametrize
122. ods 2 12 4 Class Image Include file TeF Image hxx Constructors Image const int aDim1 const int aDim2 const bool aDebug Image const Image amp alm Public methods int int int int int int int int int int int int int int Resize const int aDimi consintinimpdim2 Matrix lt double gt amp aMap ReadMapFile const string aFileName Matrix lt double gt amp aMap Import Matrix lt double gt amp aMap const Vector lt double gt aProfile Import Matrix lt double gt amp aMap const Vector lt double gt amp aProfile const Vector lt double gt amp a ImportAndNormalizeProgressive Vector lt double gt aln const double aThres1 const bool ImportAndMeanProgressive Vector lt double gt amp aln const double aThres1 const bool aCont ImportAndMeanSigmaProgressive Vector lt double gt amp aln const double aThres1 const bool ImportAndPercentileThresholdProgressive Vector lt double gt aln const double aPercent BuildEnergyHisto Matrix lt double gt amp aMap BuildFrequencyProfile Matrix lt double gt amp aMap Vector lt double gt amp aProfile BuildFrequencyProfileAndSigma Matrix lt double gt amp aMap Vector lt double gt amp aProfile Vector Thresholding Matrix lt double gt amp aMap const double aThreshold MeanThresholding Matrix lt double gt amp aMap MeanAndSigmaThresholding Matrix lt double gt amp aMap const Vector lt double gt amp aMeans const V 141 int
123. oloredNoise Noise e BuM_FFTW hxx typedef enum typedef enum FFT_Complex FFT_Estimate FFT_Real FFT_Measure Type FFT_Patient FFT_Exhaustive Wisdom typedef enum typedef enum Real FFT_Forward Imaginary FFT_Backward Complex FFT_Both Modulus Direction Output e ITF Interferometer hxx typedef enum 1 Virgo Virgo_Quartz Virgo_Quartz_Yag Hanford_2k Hanford_4k Livingston GEO TAMA Detector e Spectral_Taper hxx typedef enum None Hann Welch 147 Haming Bartlett Window e Spectral_Spectral hxx typedef enum One Sided Two Sided PSDType 148 2 16 Appendix 2 VIRGO 20 kHz sampling static const int FZERO_VIRGO 20000 VIRGO noise SIGMA_NOISE SENS_NOISE sqrt FZERO 2 static const double SIGMA_NOISE_VIRGO 4 e 21 static const double SENS_NOISE_VIRGO 4 e 23 VIRGO expected sensitivity curv parameters static const double SPEND_VIRGO 1 2e 36 static const double SMIRROR_VIRGO 3 6e 43 static const double SSHOT_VIRGO 3 5e 46 static const double FCUT_VIRGO 500 static const double SPEND_QUARTZ 4 5e 39 static const double SMIRROR_QUARTZ 3 6e 43 static const double SSHOT_QUARTZ 3 5e 46 static const double FCUT_QUARTZ 500 static const double SPEND_QUARTZ_YAG 4 5e 39 static const double SMIRROR_QUARTZ_YAG 1 2e 44 static const double SSHOT_QUARTZ_YAG 3 5e 46 static const doubl
124. onic int aSamplingFrequency double aFrequency double aCoefficient int SetElectronic double aFrequency double aSigmaProcess double aSigmaProcess double aCoefficient int SetSamplingFrequency int aSampling void SetDebug bool aValue bool GetDebug int GetProcessNoise Virgo Vector lt double gt 4 aProcessNoise double GetSigmaProcess double GetSigmaProcess double GetCoefficient double GetOmega double GetMass double GetTemperature double GetQuality 87 2 7 6 2 Class GaussianLine Include file Noise_GaussianLine hxx Constructors GaussianLine GaussianLine bool aDebug int aSamplingFreq GaussianLine string aFileName Destructors GaussianLine Public methods int Reset double aF double aQ double aMass double aT double aSigmaProcessCoef double aSigmaMeasure double aMeasureCoef int Reset int aSamplingFreq double aF double aQ double aMass double aT double aSigmaProcessCoef double aSigmaMeasure double aMeasureCoef int Reset void int GetData Vector lt double gt amp aData void SetDebug bool aValue bool GetDebug 88 2 8 BuS sources of Burst events This package includes one generic virtual C class BuS_Source from which is five other classes BuS_Gaussian BuS_DampedSine BuS_GaussianSine BuS_ZM and BuS_DFM inherit The five previous classes provide methods to fill Vector of data with Gaussian damped sine and Gaussian sine waveforms as well
125. or lt double gt amp aData int GetData double aPosition Vector lt float gt amp aData double GetValue double aX int SetSigma double aSigma int SetDistance int aDistance The constructors define the d parameter of the Gaussian function and the SNR value which is used to compute the normalization factor which takes into account the kind of noise which is chosen and the SNR value The method GetData fills the data vector with a gaussian signal which is located at the relative position given by the argument aPosition This position is a fraction of the data vector size The method GetValue returns the discrete value of the Gaussian waveform The method SetSigma changes the value of the sigma of the Gaussian waveform The method SetDistance performs no action Program example include BuS_Gaussian hxx using namespace _BuL int main int arg char argc Vector lt double gt dou size double position 5 double sigma 1 double SNR 10 Gaussian signal sigma SNR signal SetSigmaNoise 1 signal GetData position dou for int j size 2 j lt size 2 2 j cout lt lt Vector dou lt lt j lt lt lt lt dou j lt lt endl ITF antenne ITF Virgo Vector lt double gt PSD 1000000 antenne GetOneSidedPSD PSD Gaussian signal_color sigma SNR sampling PSD signal_color GetData position dou 91 for int j size 2 j lt size 2 2 j cout lt lt Vec
126. or lt double gt phase Vector lt double gt data datasize Spectral psd datasize Taper Hann 1 Spectral QneSided 1 pole simple Vector lt double gt PoleFrequency 1 Vector lt double gt PoleQ 1 Vector lt double gt ZeroFrequency 0 Vector lt double gt ZeroQ 0 file open filter dat ofstream out PoleFrequency 0 100 PoleQ 0 0 gain 1 filter new IIR ZeroFrequency ZeroQ PoleFrequency PoleQ gain sampling filter gt PrintCoefficients filter gt SetGainAtFrequency 1 1 cout lt lt Gain at 10 Hz lt lt filter gt GetGainAtFrequency 10 lt lt endl cout lt lt Global Gain lt lt filter gt GetCurrentGain lt lt endl Impulse reponse for i 0 i lt datasize i if i 0 output filter gt Filter 1 else output filter gt Filter 0 file lt lt double i sampling lt lt lt lt output lt lt endl datali output file close 44 psd PSD data out int sampling psd WriteSpectrumFile out tf dat sampling filter gt GetTransferFunction module phase 1 module WriteFile mod dat phase WriteFile phase dat delete filter return 0 45 2 5 4 IIR filter using the Butterworth approximation This class provides facilities to design low pass high pass and band pass filters using the Butterworth approximation design 9 The Butterworth filter provides the best Taylor Se
127. ort time Fourier transform based analysis and wavelet or multi resolution When the Fourier transform is used time resolution depends as usual on the delay between successive windows but also on the window shape and other operations done on the signal Frequency resolution depends on the window length plus intermediate operations For wavelet transforms especially a trous algorithm which is a shift invariant adaptation of the discrete wavelet transform time resolution is fixed by the chosen wavelet whereas frequency resolution is defined by the algorithm Most useful information can be found in 25 and 26 The classes in this library is organized by analysis step I signal transformation Transforms IT trans form data mining Image III detection algorithms TBD plus several tool classes developed for a specific use but with a larger potential Cluster and Histo 2 12 2 Class Transforms Include file TeF _Transforms hxx Constructors Transforms const int aSignalFreq const int aTailleIn const bool aDebug Public methods int int int int int int int int int int int int int int Analytic Vector lt double gt amp aData GaussianWindow const int aTransformLength const double aGaussianTO HanningWindow const int aTransformLength const double aFactor GaborParametrize const int aTransformTLength const int aTransformTStep Gabor Vector lt double gt amp aData Matrix lt double gt amp aMap Aut
128. ow double aFrequencyPassLow double aFrequencyPassHigh double aFrequencyStopHigh double aTFStopLow double aTFPassBand double aTFStopHigh Destructors Butterworth Public methods double Filter double alnput int Filter Vector lt double gt amp alnput Vector lt double gt amp a0utput int GetTransferFunction Vector lt double gt amp aModulus e Vector lt double gt amp aPhase double aDeltaFrequency void PrintCoefficients SetDebug bool aValue bool GetDebug The first constructor is used for defining a low and a high pass filter depending on the value of the stop band and pass band frequencies see Fig 2 14 The second constructor is used to define a pass band filter In all cases the filters consist to cascade second order butterworth low pass filter The order of the resulting low pass filter is computed from the arguments given in the constructor The Filter method apply on a data input the time domain filter The GetTransferFunction method computes the amplitude and the phase of the transfer function using Eq 2 5 One has to specify the frequency step to fill out the vectors The PrintCoefficients method print out the expression of all the second order filters defined and used in cascade for implemented the requested filter as follow multi pole and zero filter described in section 2 5 2 Example Pass band filter using butterworth class include lt stdio h gt include lt stdlib
129. packages including some of the external packages are managed using CMT it is necessary to install the CMT package first 1 if it is not already present BuL is making use of the external packages e FFTW v3r0p1 e DCPP v2r1 All the BuL packages can be downloaded from the LAL CVS repository or from the BuL Web site http www lal in2p3 fr recherche virgo BuL home html The architecture of this distribution is the following DataAnalysis BuL BuM v3r5 DataAnalysis BuL BuS v3r5 DataAnalysis BuL BuF v4r4 DataAnalysis BuL BuT v3r5 DataAnalysis BuL BuC vir3 DataAnalysis BuL ITF v3r4 DataAnalysis BuL LinearFilter v2r4 DataAnalysis BuL Noise v3r5 DataAnalysis BuL Spectral v2r4 DataAnalysis BuL Algorithm v2r4 DataAnalysis BuL TeF vir8 DataAnalysis BuL BuLSet v3r4 The installation of the BuL library is then easy following the instructions gt cd BuL BuLSet vXrY cmt gt cmt config gt cmt broadcast cmt config gt cmt broadcast gmake The last command will send a broadcast command gmake to all the BuL packages in order to compile them 1 3 C code style e The use of pointers at the level of the users program is source of troubles errors and memory leakage problems Moreover in most of the cases it is not needed that is why the BuL library defines and make use of objects in such a way that the user has never to manipulate pointers on objects Those objects may be simple as basic containers or much more elabor
130. paralell The method Get Window int alndex returns the size of the window corresponding to the index alndez 132 2 11 8 Class PC The class PC implements the algorithm described in 22 This a Wiener filtering using Gaussian form templates The algorithm is coded in the public method correlation This function computes the correlation function between the data vector of size N and some gaussian templates of width o which are generated and stored in memory at the configuration time The output result is returned in a Matrix of size TemplateNumber x N The filter output follows a gaussian probability distribution Note that the output for each template is normalized and then the threshold level is identical for all the templates The test application source file is test BuF_PC cxx Include file BuF_PC hxx Constructors PC string aName Public methods int Configure string aFileName int Configure Vector lt double gt amp aSigma double aThreshold int aSampling int aTemplateSize int Amplitude Vector lt double gt amp aData Matrix lt double gt amp aAmplitude int Threshold Matrix lt double gt amp aAmplitude list lt Event gt amp aEventList int Threshold Matrix lt double gt amp aAmplitude double amp aGPsS list lt Event gt amp aEventList int Detect Vector lt double gt amp aData Matrix lt double gt amp aAmplitude list lt Event gt amp aEventList int Detect Vector lt double gt amp
131. ple include BuS_ZM hxx include DataContainer_Vector hxx using namespace _BuL int main int arg char argc int i j int size 2010 double position 0 3 Vector lt double gt dou size ZM signal alb2g1 10 signal GetData position dou cout lt lt Signal SNR 10 lt lt endl for j 0 j lt 20 j cout lt lt Vector dou lt lt j lt lt lt lt dou j lt lt endl signal SetDistance 10 signal GetData position dou signal SetSNR 10 signal GetData position dou cout lt lt Signal SNR 10 lt lt endl for j 0 j lt 20 j cout lt lt Vector dou lt lt j lt lt lt lt dou j lt lt endl 98 2 8 6 Class DFM The Dimmelmeier amp Font amp Mueller catalogue gathers a collection of waveforms produced by the relativist simulation of the collapse of a rotating stellar core to a neutron star assuming an axisymetry The GW field amplitude emitted during the core collapse is derived from the quadrupole moment By symetry 2D and axisymetry there is only one non vanishing term A AZ The GW amplitude is then given by 1 15 Az 24 gt 205 2 2 hy 3 sin Oe 2 26 hy 0 The DFM catalogue contain 25 waveforms obtained for different value of the following parameters e A Angular momentum eB f me ratio between the rotational and the potential energies e G IT reduced adiabatic indice The Table 2 8 6 gives the initial values
132. r according to the following procedure The sort in the increasing start time order the list of event combine any ordered couple i j of events if Ti o0 gt To see The 2 events are then combined to define a new event k Let s asume that 7 is the event index whose amplitude at maximum is the highest then we have gt Tejon Totari 7 ee e me max max es T naz 7 ak o Template Numero i Template Number 1 Micro event list contains the 2 events remove the events i and j from the list of events add the event k in the list of events continue to cluster the modified list of events ClusterCorrelatorLoose is identical to the ClusterCorrelatorTight algorithm except that one defines and uses internally the end time of a macro event as the maximum of the 2 micro events end time and the one of the maximum amplitude template This enlarges the clusterization The PushBackMicroEvent method does a push back of an event into the list of micro events 126 The PushFrontMicroEvent method does a push front of an event into the list of micro events The ClearMicroEvent method clears the list of micro events The GetMicroEvent returns a copy of the list of micro events o clusterization g i O no clusterization o Figure 2 22 Clusterization of 2 micro events top and example of 2 micro events which do not fulfill the criteria to be clusterized bottom
133. r and rectangle coefficients The SetLength double aValue method sets the Length to aValue The SetWidth double aValue method sets the Width to aValue The SetTiltAngle double aValue method sets the Tilt Angle to aValue The GetLength method returns the Length The GetWidth method returns the Width The GetTiltAngle method returns the value of the tilt angle between the x axis and a large side The Characteristic TilerPoint amp aPoint allows to test if aPoint belongs or not to the rectangle The Area returns the rectangle area 114 2 10 6 Class Tiler This class a virtual one which can not be instancied It will be inherited by all tiling procedures Include file Tiler hxx Constructors Tiler int aDim 2 Destructors Tiler Public methods virtual int TileQ 0 int GetDim void SetPrintPeriod int aPeriod int GetPrintPeriod void SetDecreaseSurfaceFactor double aValue double GetDecreaseSurfaceFactor void SetDebugLevel int aValue int GetDebugLevel void SetBoundMin int alndex double aValue double GetBoundMin int alndex void SetBoundMax int alndex double aValue double GetBoundMax int alndex void SetCharacteristicFunction bool aFunction TilerPoint amp aPoint void SetMarginCharacteristicFunction bool aFunction TilerPoint amp aPoint void SetAreaFunction double aFunction void SetSpaceParamArea double aValue double GetSpaceParamArea double
134. ries approximation to the ideal lowpass filter response at analog frequencies Butterworth low pass filters are defined by the property that the magnitude response is maximally flat in the pass band For a N order low pass filter this means that the firts 2N 1 derivatives of the magnitude squared function are zero at Q 0 Another property is that the magnitude squared function is monotonic in the pass band and the stop band The magnitude squared function for a continuous time Butterworth low pass filter is of the form 1 H GQ 27 OOP a 2 7 Using s jQ Eq 2 7 writes 1 H s H s r 2 8 c s H s 1 Y 2 8 The roots of the denominator polynomial are therefore located at values of s Sp N e0TRMORD k 0 1 2N 1 2 9 There are hence 2N poles equally spaced in angle on a circle of radius Qe in the s plane A pole never falls on the imaginary axis On the contrary if N is odd poles exists at Q and Q Finally note that poles always appear in complex conjugated pairs and that a poles at sx corresponds a pole at s That implies that to build H s one just needs to consider the N poles located on the left half plane part of the s plane see Fig 2 7 That insures the filter to be stable and causal Figure 2 7 Location of the 8 poles of a Butterworth low pass filter of order N 4 The continuous time transfer function writes e N is even H s gt 2 10 e N is odd
135. rm one dimensional Forward and Backward Discrete Fourier Transform DFT of real and complex discrete data sets The Forward DFT of an one dimensional arry Xy of size N is given by pi 2Hijk j 0 The Backward DFT output is normalized such as the consecutive actions on a vector DFT Forward gt DFT Backward provides the identical input vector 2 2 1 2 FFTW description The FFTW Fourier Transform provides an adaptative implementation of the Cooley and Tukey Fast Fourier Transform To perform the FFT FFTW creates a plan which is a data structure containing all the information that FFTW needs to compute the FT The plan provides the best combination of codelets which computes FT of a given size In other words the size N of the FT is factorized such as the resulting recursive computation of smaller F T gives the fastest answer This factorization is called WISDOM and can be saved on disk to be reused later on This plan is created at runtime and is adapted to the hardware architecture The definition of a plan depends basicaly on the size of the input data vector and on the FT direction This particulary means that a plan can be must be reused as many times as identical FT have to be computed Different kinds of plan determination wisdom can be chosen e ESTIMATE the factorization is guessed without measure The result is immediate but not optimal e MEASURE the executation time of many plans is measured with the pragmatic hypothesi
136. rns the value of the ellipse small axis The GetTiltAngle method returns the value of the tilt angle between the x axis and the long axis The Characteristic TilerPoint amp aPoint allows to test if aPoint belongs or not to the ellipse The Area returns the ellipse area 113 2 10 5 Class TilerRectangle This class inherits from TilerPattern It defines the rectangular pattern As TilerPattern it is only for internal use of the Tiler package No use should be done by standard users Include file TilerRectangle hxx Constructors TilerRectangle double aLength 1 double aWidth 1 double aTilt 0 TilerRectangle const TilerRectangle amp aRectangle Destructors TilerRectangle Public methods void SetLength double aValue void SetWidth double aValue void SetTilt double aValue double GetLength double GetWidth double GetTiltAngle bool Characteristic TilerPoint aPoint double Area A rectangle is defined by its Center its Length its Width and its Tilt Angle defined as the agnle between the x axis and a large side The constructor TilerRectangle double aLength 1 double aWidth 1 double aTilt 0 creates an rectangle centered on the origin with Length aX 2Coef f Width aY 2Coef f and Tilt aTilt By default an unit square is created The constructor TilerRectangle const TilerRectangle amp aRectangle creates an rectangle which is an exact copy of aRectangle cente
137. rray is used to fill the starting point coordinates The GetStartingCenter method returns the starting point of the tiling procedure 116 2 10 7 Class TilerWithEllipses Include file TilerWithEllipses hxx Constructors TilerWithEllipses Destructors TilerWithEllipses Public methods void SetX2Function double aFunc TilerPoint amp aPoint void SetY2Function double aFunc TilerPoint amp aPoint void SetXYFunction double aFunc TilerPoint amp aPoint void SetSafetyFactor double aValue double GetSafetyFactor void SetAngleUffset double aValue double GetAngleUffset void SetTurnDirection double aValue double GetTurnDirection void SetMinimalCoveredSurface double aValue double GetMinimalCoveredSurface void SetNbLoopsComputeCoveredSurface int aValue int GetNbLoopsComputeCoveredSurface int Tile double TotalEllipsesArea double EstimatedCoveredArea void WriteKumac std string aFileName int aPeriod The constructor TilerWithEllipses creates the object No parameter is needed The SetX2Function double aFunc TilerPoint amp aPoint method allows to define the function which returns the value of the X2Coeff for the ellipse centered on aPoint The SetY2Function double aFunc TilerPoint amp aPoint method allows to define the function which returns the value of the Y2Coeff for the ellipse centered on aPoint The SetX YFunction double aFunc
138. ructors Goertzel bool aDebug Destructors Goertzel Public methods int Setup double aSamplingFrequency double amp aLineFrequency int amp aChunkSize double Filter Vector lt double gt amp aData int Filter Vector lt double gt amp aData Vector lt double gt amp aLine int Filter Vector lt double gt amp aData Vector lt double gt amp aLine Vector lt double gt amp aWindow void Reset The method Setup requires the values for fo f and N aChunkSize Once called the filter memory is initilized to 0 Filtering either produces one DFT coefficient for input Vector of size larger than N or a Vector of values for the coefficient computed on non overlapping chunks of length N It is also possible to apply a window on the input data The use of these additional operations is not very clear since thefilter deals with one frequency only at a resolution fixed by the size of the non null 59 part of the data The Reset method re initializes the filter memory to 0 Program exemple include LinearFilter_Goertzel hxx include Noise_Gaussian_MT hxx include Noise_GaussianAR hxx include Spectral_Spectral hxx include Spectral_Taper hxx include DCPP_Matrix hxx using namespace _BuL AAA AAA AAA AAA AAA AAA AAA AAA AAA AAA LT int main int argc char argv variables Goertzel filtre false double samplingFrequency 20000 int dataSize int pow 2 0 20 Vector lt dou
139. s nominal fo 16384 Hz sampling fre quency The other constructor allow to choose a different sampling frequency if needed The methods are described in the section of the ITF class Program exemple include ITF_GE0O hxx using namespace _BuL int main int arg char argc int i j int size 10000 Vector lt double gt psd_sdt size GEO std 10000 std GetOneSidedASD size psd_sdt for i 0 i lt size 1000 i cout lt lt vector ASD lt lt i lt lt lt lt psd_sdt i lt lt endl return Q 70 2 6 6 Class TAMA Include file TTF_TAMA hxx Constructors TAMA TAMA int aSamplingFrequency Public methods int GetTwoSidedPSD int aSize Vector lt double gt amp aData int GetOneSidedPSD int aSize Vector lt double gt amp aData int GetTwoSidedASD int aSize Vector lt double gt amp aData int GetOneSidedASD int aSize Vector lt double gt amp aData void SetSamplingFrequency int aValue int GetSamplingFrequency double GetLatitude double GetLongitude double Get Azimuthal double GetArmAngle void SetLatitude double aValue void SetLongitude double aValue void SetAzimuthal double aValue void SetArmAngle double aValue void SetDebug bool aValue bool GetDebug The TAMA constructor consider the ITF TAMA detector with its nominal fo 16384 Hz sampling frequency The other constructor allow to choose a different sampling frequency if needed The methods are descr
140. s returned in a vector filled up to bin N n 1 The end of the vector from bin N n up to bin N 1 is filled with zeros The filter output follows a Chi Square probability distribution The filter output is given by X X 2a Xa Xo lz 3 n 1 a cov Xa K Na aed where X and X are the normalised slope and the offset of the linear regression calculated in the moving window of size n ALF The test application source file is test BuF_ALF cxx Include file BuF_ALF hxx Constructors ALF string aName Public methods int Configure string aFileName int Configure Vector lt int gt amp aWindowSize Vector lt double gt amp aThreshold int aSampling int Amplitude Vector lt double gt amp aData Matrix lt double gt amp aAmplitude int Threshold Matrix lt double gt amp aAmplitude list lt Event gt amp aEventList int Threshold Matrix lt double gt amp aAmplitude double amp aGPS list lt Event gt amp aEventList int Detect Vector lt double gt amp aData Matrix lt double gt amp aAmplitude list lt Event gt amp aEventList int Detect Vector lt double gt amp aData double amp aGPS list lt Event gt amp aEventList int End int SetThreshold Vector lt float gt aThreshold int GetThreshold Vector lt float gt amp aThreshold int GetNWindow int GetWindow int alndex void SetDebug bool avalue bool GetDebug int Regression Vector lt double gt amp aData i
141. s that if an optimal plan for a size N is known this plan is still optimal when size N is used for a larger transform The drawback of this feature is that it can take a long time from seconds to minutes depending on the machine and the precision of the clock timer to create the plan hence it is recommended only when one repeats FT request e PATIENT improved algorithm compared to MEASURE especially adapted for long array but longer e EXHAUSTIVE improved algorithm compared to PATIENT This is the most optimal wisdom but the longest minutes and such should be used when a plan is used many times 1FFTW 3 0 1 since release BuM v4 2according to the execution time and not the number of floating point operations 3which are small piece of code generated at compile time A wisdom can be saved on disk and then reused for any size which it is applicable as long as the planner option choice is not more demanding than this with which the wisdom has been cerated For example a MEASURE wisdom can be used only for ESTIMATE and MEASURE but not for PATIENT and EX HAUSTIVE Finally FFTW implements a DFT for complex data as well as a DFT for real data set The latter gains roughly a factor 2 for the computation time and storage 2 2 1 3 Class FFT description The class FFT embeds the FFTW C library providing its main functionalities Five different constructors have been designed in order to fulfill different working cases As minimal informatio
142. set the two sided PSD which is used for performing the match filtering Note that the template is then re normalized using the new PSD The method SetTemplate set the template to the content given as argument Note that the tem plate will be normalized using the current PSD given at object creation or using the SetPSD method The vector containing the template must have a size equal or smaller than the size given in the constructor The method Filter performs the match filtering between the input data vector and the template The result is stored in the output vector The method GetNormalizationFactor returns the factor used for the template normalization Program example include math h include LinearFilter_MatchFilter hxx using namespace _BuL int main int narg char argc int i double x double sigma 10 int size 10000 int sampling 20000 Vector lt double gt patron size Vector lt double gt data size Vector lt double gt output size Vector lt double gt psd_twosided 2 size 97 double factor double norm sqrt sqrt M_PI sigma for i 0 i lt size i x i size 2 patron i exp x x 2 sigma sigma norm psd_twosided ReadFile data psd dat MatchFilter filter sampling psd_twosided patron filter SetDebug true factor filter GetNormalizationFactor for i 0 i lt size i patron i patron i factor
143. ssianSine double aSigma double aFrequency double aSNR int aSampling Vector lt double gt amp aPSD GaussianSine double aSigma double aFrequency double aSNR ITF Detector aName int GetData double aPosition Vector lt double gt amp aData int GetData double aPosition Vector lt float gt amp aData double GetValue double aX int SetSigma double aSigma int SetFrequency double aFrequency int SetDistance int aDistance The constructors define the o and f parameters of the GaussianSine function and the SNR value which is used to compute the normalization factor which takes into account the kind of noise which is chosen and the SNR value The method GetData fills the data vector with a GaussianSine signal which is located at the relative position given by the argument aPosition This position is a fraction of the data vector size The method GetValue returns the discrete value of the GaussianSine waveform The method SetSigma changes the value of the o of the GaussianSine waveform The method SetFrequency changes the value of the f parameter of the GaussianSine waveform The method SetDistance performs no action Program example include BuS_GaussianSine hxx using namespace _BuL int main int arg char argc double position double sigma 5 0 1 double frequency 300 double SNR 10 GaussianSine signal sigma frequency SNR signal SetSigmaNoise 1
144. st double aScale Pix GetPix const int aNum Point GetGeoBary void Point GetWeightedBary void 139 Point GetPeakPoint void int GetSize void double GetWeight void double GetMean void double GetVariance void double GetDimiSpread void double GetDim2Spread void double GetDimiSpread double amp aMin double amp aMax double GetDim2Spread double amp aMin double amp aMax double GetM11 void double GetM02 void double GetM20 void double GetEqRadius void int GetPeakNb void double GetPeakPower void bool GetCloseToBorder void int GetImageNum void void AppendFile const string aFileName const int aRefNum void SetDebug const bool aValue bool GetDebug The cluster class provides a list of Pix along with a list of methods for parameter extraction Some methods are also used within the clustering algorithm to create expand and merge clusters as required The Pix structure holds two pixel integer coordinates fAxis1 fAxis2 and a value fVal double It is used to store pixels of a cluster Some operations on Pix have been implemented PixDistance that re turns the Euclidean distance that checks identity and lt checking strict inferiority of both coordinates The Point structure holds two map double coordinates fX1 fX2 and a value fVal double It is used to store the geometrical barycentre of a cluster One operation is available PointDistance similar to PirDistanc
145. structor consider the ITF Livingston detector with its nominal fo 16384 Hz sam pling frequency The other constructor allow to choose a different sampling frequency if needed Note that for the moment the spectrum power is identical for the 3 LIGO detectors The methods are described in the section of the ITF class Program exemple include ITF_Livingston hxx using namespace _BuL int main int arg char argc int i j int size 10000 Vector lt double gt psd_sdt size Livingston std 10000 std GetOneSidedASD size psd_sdt for i 0 i lt size 1000 i cout lt lt vector ASD lt lt i lt lt lt lt psd_sdt i lt lt endl return 0 69 2 6 5 Class GEO Include file TTF_GEO hxx Constructors GEO GEO int aSamplingFrequency Public methods int GetTwoSidedPSD int aSize Vector lt double gt amp aData int GetOneSidedPSD int aSize Vector lt double gt amp aData int GetTwoSidedASD int aSize Vector lt double gt amp aData int GetOneSidedASD int aSize Vector lt double gt amp aData void SetSamplingFrequency int aValue int GetSamplingFrequency double GetLatitude double GetLongitude double Get Azimuthal double GetArmAngle void SetLatitude double aValue void SetLongitude double aValue void SetAzimuthal double aValue void SetArmAngle double aValue void SetDebug bool aValue bool GetDebug The GEO constructor consider the ITF GEO detector with it
146. t endl lt lt endl Hanford H2k ITF Hanford_2k 10000 H2k GetOneSidedASD size psd_2k for i cout double double double double cout lt lt cout lt lt cout lt lt cout lt lt return O i lt size 1000 i lt lt vector ASD lt lt i lt lt lt lt psd_2k i lt lt endl longitude H2k GetLongitude latitude H2k GetLatitude azimuthal H2k GetAzimuthal arm_angle H2k GetArmAngle Hanford Longitude lt lt longitude lt lt endl Hanford Latitude lt lt latitude lt lt endl Hanford azimuthal angle lt lt azimuthal lt lt endl Hanford arm angle lt lt arm_angle lt lt endl 0 68 2 6 4 Class Livingston Include file TTF_Livingston hxx Constructors Livingston Livingston int aSamplingFrequency Public methods int GetTwoSidedPSD int aSize Vector lt double gt amp aData int GetOneSidedPSD int aSize Vector lt double gt amp aData int GetTwoSidedASD int aSize Vector lt double gt amp aData int GetOneSidedASD int aSize Vector lt double gt amp aData void SetSamplingFrequency int aValue int GetSamplingFrequency double GetLatitude double GetLongitude double GetAzimuthal double GetArmAngle void SetLatitude double aValue void SetLongitude double aValue void SetAzimuthal double aValue void SetArmAngle double aValue void SetDebug bool aValue bool GetDebug The Livingston con
147. t GetZerosNumber The ComputeARMA method computes the AR and MA coefficients from the PSD given as argument using the procedure described above The aPSDScaleFactor argument is the scale factor which has been applied to the PSD in order to avoid to manipulate very small number The aMinimalFrequency argument is the frequency below which the PSD is taken as constant and equal to the value at aMinimalFrequency The ComputeAR method computes the AR and MA coefficients as the previous method but taking into account the official PSD which corresponds to the given interferometer The aMinimalFrequency argument is the frequency below which the PSD is taken as constant and equal to the value at aMini malFrequency The SaveARMA method saves in the file whose name is given as argument the ARMA filter coeffi cients according to the format which is described in section 2 7 3 2 The GetAR and GetMA methods return the values of the coefficients of the AR and MA part of the ARMA filter The SetDeltaFrequency and GetDeltaFrequency methods allow to change and to get the value of the frequency step which is used when computing the autocorrelation function The SetPolesNumber and GetPolesNumber methods allow to change and get the number of poles of the AR filter The SetZerosNumber and GetZerosNumber methods allow to change and get the number of poles of the AR filter The SetSamplingFrequency
148. t double gt amp GetWaveForm Matrix lt complex lt double gt gt amp GetFourier WaveForm When one creates a TGaussian object the number of templates are computed This number depends on the 3 parameters e Omin minimal value of the templates gaussian width e Omaz Maximal value of the templates gaussian width e e maximal value of the SNR loss The minimal match MM is defined MM 1 e The number of templates is given by Tmin dog lo maz y one 4 4 where y and y_ are given by VI MM2 V1 MM 1 MM MM The different values of gaussian template width are determined using the recursive relation Oi 1 Sa 05 00 Omi 1 1 1 y i 0 min 102 The method Generate WaveForm fills a template matrix with values of Gaussian given in 2 27 The signal is centered at the middle of the data vector The method GenerateFourierWaveForm fills a template matrix with Fourier Transform of the time domain Gaussian templates One just applies a direct FFT on th etime domain wave form The methods GetWaveForm and GetFourierWaveForm return the pointer to the vector containing the template generated in time or frequency space One should mention that all templates are stored in vector of equal size given by the largest template Omaz The size of the vector containing the templates in time domain is the power of 2 just greater than 8 X Omar X Sampling where Sampling is supposed to be 20 kHz The size
149. tData Vector lt double gt amp aData int GetData Vector lt float gt amp aData double Get Value int SetSamplingFrequency int aValue int GetSamplingFrequency int SetGaussianSigma double aSigma The class Gaussian MT is making use of the Random_MT class to generate random numbers The first constructor initializes the sigma of the white and Gaussian noise to the official Virgo value The second constructor allows to specify the sigma value The GetData methods return a vector of the given size of values generated according to a white noise random distribution Note that the memory allocation of the vector has to be done by the user This method makes use of the Random_MT GetGaussianData and hence there is no change of Seed number before generating data It is up to the user to do it by calling SetSeed or SetSeedLocalTime methods The GetValue method returns a scalar value generated according to a white noise random distribution The conditions are identical to the GetData method The SetSamplingFrequency method changes the sampling frequency private variable and also the sigma value of the noise according to the formula given above The sampling frequency is used only in the computation of the sigma when one uses the constructor I The GetSamplingFrequency method gets the value of the sampling frequency used On the contrary the SetGaussianSigma method changes only the value of the sigma of the no
150. ta vector size The method GetValue returns the discrete value of the DampedSine waveform The method SetTau changes the value of the 7 of the DampedSine waveform The method SetFrequency changes the value of the f parameter of the DampedSine waveform 93 The method SetDistance performs no action The method SetDamped allows to switch from one quadrature to anothre at any time The change is taken into account at the next call of GetData Program example include BuS_DampedSine hxx using namespace _BuL int main int arg char argc Vector lt double gt dou size double position 25 double tau 0 005 double frequency 300 double SNR 10 int damped 1 DampedSine signal tau frequency damped SNR signal GetData position dou for int j size 4 16 j lt size 4 19 j cout lt lt Vector dou lt lt j lt lt lt lt dou j lt lt endl DampedSine signal_itf tau frequency 10 ITF Virgo signal_itf SetSNR 100 signal_itf GetData position dou for int j size 2 j lt size 2 2 j cout lt lt Vector dou lt lt j lt lt lt lt dou j lt lt endl return 0 94 2 8 4 Class GaussianSine The class GaussianSine provides facilities to generate a waveform Include file Constructors Public methods 42 f t ae 27 sin 27 ft BuS_GaussianSine hxx GaussianSine double aSigma double aFrequency double aSNR Gau
151. the wisdoms corresponding to forward and backward FT plan into 2 files The file names are Forward fileName and Backward fileName Note that when creating a FFT object giving the fileName corresponding to the wisdoms the construc tors IV amp V follow the file name syntax Forward fileName and Backward_fileName The method GetSize gives the size of the plan hence the size of the data vector which are processed by the Fourier transform Program example Complex FT int size 1024 Vector lt double gt input size Vector lt complex lt double gt gt output size Vector lt complex lt double gt gt end size FFT A size FFT FFT_Complex A Forward input output A Backward output end 5the modulus of a complex number X is defined as VX X 11 Program example Real FT int size 1024 Vector lt double gt input size Vector lt complex lt double gt gt output size Vector lt double gt end size FFT B size FFT FFT_Real B Forward input output B Backward output end 12 2 2 2 Class Correl This class provides several methods to perform special mathematical processing such as the correlation the autocorrelation and the convolution of discrete data sets All the methods are based on the Fourier Transform and hence make use of the FFT class The main interest of the Correl object is due to the fact that the FFTW plans are determined only once when the o
152. ther classes are derived from this base class 2 4 1 1 Class MovingSlice Let s assume one wants to compute a quantity with a sub set of data corresponding to a window size w and that one wants to repeat it moving the window by step of A bins parameter definition constraint input data size N window size W w lt N slice step A O lt A lt N output data size x x is an integer memory size w A In the MovingSlice class one defines which part of the data chunk is used to compute the element of the output vector of the algorithm As shown in Figure 2 2 e element 0 input data the w A last data of the previous data chunk the A first data of the input data chunk e element 1 input data the w 2A last data of the previous data chunk the 2 A first data of the input data chunk e element x 1 input data the w last data of the data chunk initial memory null memory memory rr aes A ee OUTPUT Figure 2 2 Slicing mechanism of the MovingSlice algorithm 32 The content of the class MovingSlice is given for the sake of completion but one can not instance such an object directly It is used by other classes such as the Burst filters There is only one public method This content corresponds to the minimal content of any Algorithm class Include file Algorithm MovingSlice hxx Constructors MovingSlice int alnputSize int aWindowSize int aMovingDelta MovingSlice int alnputSize int
153. time frequency w associated to the z variable z e is Due 2i A ae Gtan 5 s 0 1 A 2 hence o 0 and Q lt tan The bilinear transform is a method of squashing the infinite straigh analog frequency axis so that it becomes finite But to avoid squashing the filter s desired frequency response too much the bilinear transform squashes the far ends of the frequency axis the most leaving the middle part relativelty unsquashed as shown in Fig 2 3 2 2 gr tan 02 itil RD 43 03 Ura Ue Figure 2 3 Frequency warping induced by the bilinear transfrom of a continuous time pass band filter H jQ into a discrete time band pass filter H e Extracted from 10 This effect related to the non linear property of the bilinear transform changes the shape of the desired filter frequency response especially in the transition bands region It is acceptable when designing filters 38 containing piecewise constant frequency magnitude characteristics such as High pass Low Pass and Band Pass filters 2 5 2 How to design a linear filter with class IIR using the zeros and poles design The IIR class prevides an implementation of LTI digital filter providing the poles and zeros in continuous time domain It allows to define filters composed of numerous pass bands and attenuation bands whose caracteristics are constant One has to specify the zeros and the poles of the filter following the convention given in section
154. tion 2 15 10 Type FFT_Complex Output Real FFT_Real Imaginary Complex Modulus Wisdom FFT_Estimate Direction Forward FFT Measure Backward FFT Patient Both FFT Exhaustive The Forward method performs a forward FFT of a Vector of real float or double or complex float or double values The output is stored in a complex Vector or in a real Vector half complex format depending on the choice of the user The Backward method will perform an backward FFT of a Vector of complex values or of input Vcetor of half complex values The output is stored in a complex Vector When dealing with real Vectors it is better to use the real FFT defined in FFTW gain of a factor 2 in computional time as shown in the following test programs example The method SetCheckVectorSize allows to enable disable the memory size of the input and output vectors check which is performed before any FT The WriteFile methods prints out in a file the content of the output vector of a FT associated to the value of the corresponding frequency In the case of a complex Vector output we can choose to print out the full complex representation or the real imaginary or modulus part of the vector The ExportLastWisdomFile char fileName method allows after creating a FFT object to save the last Wisdom into a file Note that will correspond to the backward FT plan The SaveWisdomFile char fileName method allows after creating a FFT object to save
155. tor dou lt lt j lt lt lt lt dou j lt lt endl return 0 92 2 8 3 Class DampedSine The class DampedSine provides facilities to generate a waveform f t ae sin 2r ft and f t 0 e Fcos 27 ft Include file BuS_DampedSine hxx Constructors DampedSine double aTau double aFrequency int aDampedType double aSNR DampedSine double aTau double aFrequency int aDampedType double aSNR int aSampling Vector lt double gt amp aPSD DampedSine double aTau double aFrequency int aDampedType double aSNR ITF Detector aName Destructors DampedSine Public methods int GetData double aPosition Vector lt double gt aData int GetData double aPosition Vector lt float gt aData double GetValue double aX int SetTau double aTau int SetFrequency double aFrequency int SetDistance int aDistance int SetDamped int aDampedType The constructors define the r and f parameters of the DampedSine function and the SNR value which is used to compute the normalization factor which takes into account the kind of noise which is chosen and the SNR value The aDamped argument defines the chosen quadrature as follow e 0 cosine e 1 sine Remind that one has to give the one sided PSD in case of a colored noise The method GetData fills the data vector with a DampedSine signal which is located at the relative position given by the argument aPosition This position is a fraction of the da
156. tor is not filled up The method SetDFM changes the waveform by a new one The method SetDistance initializes the distance and the normalization factor is recomputed Program example include BuS_DFM hxx include DataContainer_Vector hxx using namespace _BuL int main int arg char argc int i j int size 2010 double position 0 3 Vector lt double gt dou size DFM signal alb2gi 10 signal GetData position dou cout lt lt Signal SNR 10 lt lt endl for j 0 j lt 20 j cout lt lt Vector dou lt lt j lt lt lt lt dou j lt lt endl signal SetDistance 10 signal GetData position dou signal SetSNR 10 signal GetData position dou cout lt lt Signal SNR 10 lt lt endl for j 0 j lt 20 j cout lt lt Vector dou lt lt j lt lt lt lt dou j lt lt endl 100 2 8 7 Class Signal Include file BuS_Signal hxx Constructors Signal Signal Vector lt double gt amp aPSD int aSampling Signal double aSigmaNoise int aSampling Signal ITF amp alTF Destructors Signal Public methods int SetSignal Vector lt double gt amp aData int SetGaussianSignal double aSigma int SetGaussianSineSignal double aSigma double aFrequency int SetGaussianCosineSignal double aSigma double aFrequency int SetDampedSineSignal double aTau double aFrequency int SetDampedCosineSignal double aTau double aFrequency
157. tring amp aKey int amp alndex string amp aValue void WriteCard void SetDebug bool aValue bool GetDebug The ReadCard opens the file reads and stored all keyword associated line in a list of strings The method returns 0 except when a problem occured during the file open phase return code 1 The GetParam searchs for a value corresponding to the given position in the line associated to the given keyword It is asumed that a keyword is unique In the following example one reads the example data card using the GetParam methods Example PARAMETERS DEFINITION ALGO Example CLUSTER Base FILTERMEMORY No NWINDOW 10 WINDOWS 10 20 30 40 50 100 150 200 250 300 THRESHOLD 31 26 30 84 30 66 30 70 30 42 29 98 29 58 28 68 28 54 28 20 6 6 k Program example include BuM_DataCard hxx using namespace _BuL int main int narg char argc int rc 22 int nWindow double seuil string name double test DataCard dc test ConfigFile txt dc SetDebug true dc ReadCard rc dc GetParam NWINDOW 0 nWindow cout lt lt nWindow lt lt endl rc dc GetParam WINDOWS 6 nWindow cout lt lt nWindow lt lt endl rc dc GetParam WINDOWS 10 nWindow cout lt lt nWindow lt lt endl rc dc GetParam WINDOWS 9 nWindow cout lt lt nWindow lt lt endl rc dc GetParam THRESHOLD 8 seuil cout lt lt seuil lt lt endl
158. ual to 2 x aSize the input vector must contain the aSize values followed by aSize 0 although the output vector size must also be equal to aSize In the case of no zero padding use aSize is the dimension of the input and output vectors The correlation function of two functions g t and h t is defined as a he g r f g t 7 h t dt The correlation theorem says heg f MP The convolution function of two functions g t and h t is defined as o e i jemnou The convolution theorem says hog f h f a f The practical implementation of the convolution and correlation functions are making use of the previ ously given properties The output vector is correctly ordering taking into account the fact that the FFT algorithm provides a vector in a wrap around order In the case of zero padded vectors the output vector is correctly formatted spoiled data part removed Program example include math h include BuM_Correl hxx using namespace _BuL int main int narg char argc int i double x double sigma 1 int size 32 Vector lt double gt a 2 size Vector lt double gt b 2 size Vector lt double gt c size Vector lt complex lt double gt gt aa 2 size Vector lt complex lt double gt gt bb 2 size Vector lt complex lt double gt gt cc size for i 0 i lt size i x i size 2 5 ali exp x x 2 sigma sigma bli a i aa i
159. ze list lt Event gt eventList Example created Example algo Example algo setDebug true Generator created double sigmaNoise 1 128 Gaussian_MT Noise sigmaNoise Noise GetData h Source double sigma 0 001 double SNR 20 double position 5 Gaussian Signal sigma SNR Signal GetData position g for i 0 i lt size i h i h i g i algo Configure test Example txt algo Amplitude h Amplitude int eventNumber algo Threshold Amplitude eventList if eventNumber gt 0 cout lt lt Micro events found lt lt endl cout lt lt eventList lt lt endl clusterization list lt Event gt macroEvent ClusterBase eventList cout lt lt Macro events found lt lt endl cout lt lt macroEvent lt lt endl Alternative method algo Detect h Amplitude eventList algo End return 0 Here follows of an example of data card used for the configuration ALF PARAMETERS DEFINITION ALGO MF CLUSTER Base FILTERMEMORY No SAMPLING 20000 NWINDOW 10 WINDOWS 10 15 20 25 30 40 50 70 THRESHOLD 5 33 5 33 5 33 5 33 5 33 5 33 5 33 5 33 OOOOHOOO 129 150 200 5 33 5 33 2 11 5 Class ALF The class ALF implements the algorithms described in 23 The algorithm is coded in the public method Regression This function computes the linear regression of N points over a moving window of size n The output result i

BuL Virgo Gravitational Wave Searches Library for Burst

Contents

Download Pdf Manuals

Related Search

Related Contents