LALAPPS_COMPUTEFSTATISTICBINARY WALKTHROUGH

A brief introduction

• The lalapps_ComputeFStatisticBinary code has been written to search for continuous GW's emmitted from sources in binary elliptical orbits.

• More specifically, the first use of the code will be to search for GW's from Sco X-1, an accreting neutron star, in an approximately circular orbit about its lower mass companion star.

• The spin frequency (and therefore GW frequency) has been constrained to a fairly wide bandwidth/bandwidths (~40 Hz) and we have therefore decided to adopt a frequency domain search strategy.
  

• The added complication to this analysis (as opposed to the isolated pulsar search) is the introduction of source motion relative to our inertial reference point , the solar system barycenter (SSB).

• A GW signal from a source in a binary system, upon interaction with a detector, will exibit Doppler shifts due to the (exactly known) motion of the detector relative to the SSB, in addition to the (not exactly known in the case of Sco X-1) orbital motion of the source.

• In addition to the Doppler modulation of the signal there will also be an amplitude modulation effect dependent upon the (known) directional response of an interferometric detector coupled with the (unknown) polarisation angle of the GW and inclination angle of the source. 
  

• We use the search code to calculate the F statistic, a quantity that maximises over the polarisation and inclination angles, allowing us to ignore them, and not include them explicitly as search parameters. JKS98

• What we do search over is a parameter space defined on the unknown orbital parameters of the system + the GW frequency. 
  

Hopefully it all seems familiar ?

From the very start of the work done on searching for neutron stars in binary sytems, and in particular Sco X-1, we have adopted the frequency domain approach. We have also deliberately tried to use (where possible) the same codes and search strategies as used in the search for GW's from isolated pulsars. Some key points and major differences between the codes are :-

1. The latest version of the binary search code is essentially identical to the isolated search code.  In fact, with the binary input flag set to false, the code is identical to the isolated pulsar search code lalapps_ComputeFStatistic.

2. We ignore the inbuilt sky position template functions. We deal at present with a single known sky position.

3. Instead of looping over an bank of sky positions we loop over a bank of ready-made orbital templates.

4. At the heart of the code we still deal with LALDemod.c but the input parameters to this function are calculated using ComputeSkyBinary.c and not ComputeSky.c.

Walkthrough

Link to the source code in the LALApps CVS

In this section we describe the code on a function call by function call basis.  Important parameter choices are mentioned as well as differences between this strategy and the isolated pulsar search strategy.  Further technical details describing each function will be discussed in future code review meetings.

• initUserVars()
• Description : Initialises the user input variables and sets them to their default values
• Key Points
  
• A binary flag is included to specify wether we do a binary search or an isolated search.
• A binary template bank file is now a required user input if the biniary flag is set.
• We have also added an option to input the running median window size (explained later).
• We have also added an option to input a variable (dopplermax) that defines how large a change in frequency can be expected from a particular source (For binaries this number is typically orders of magnitude larger than that for isolated pulsars and is dependent upon many variables). 
  

• LALUserVarReadAllInput()
• Description : Reads in user defined input variables from a configuration file or the command line.
• Key Point : If the binary flag is set then unrequired variables are read in but from this point on are ignored.

• SetGlobalVariables()
• Description : Takes the user defined input variables and generates quantities more useful to the workings of the code.
• Checks for the existence of the data and extracts information regarding the format of the data ie. SFT time baseline, frequency band, the start and end times of the data.
• Initialises data structures relating to the correct interferometer.
• Calculates the indexes of starting and final frequency bins  based on the frequency band we wish to search.
• Sets a default frequency resolution if not defined by the user.
• Allocates memory for the output F statistic.
  
• Key Points
• We have defined a different default frequency resolution (=1/5Tobs) as opposed to  ((=1/2Tobs).
• If the binary flag is set then we do not do initialise and check the sky grid parameters.
  

• AllocateMem()
• Description : Allocates memory for the majority of appropriate data structures.
  
• Key Points
• For simplicity at present we not allocate memory for spin down related structures if the binary flag is set.
• ReadSFTData()
• Description : Reads in all the data found in the directory specified by the user input that has SFT base name also specified in the user input.
• Key Point
• There has been no change to this function but we should mention at this point that we use 60 second SFT's and not the 1800 second SFT's used in the isolated search.  This is because we need to be sure that the Doppler shifted GW frequency of Sco  X-1 remains within a single frequency bin for the duration of each SFT.
• NormaliseSFTDataRngMdn()
  
• Description : This function takes each SFT, estimates the running median across the frequency band, then normalises by this running median.
  
• Key Points : Due to our shorter SFT time baseline our initial SFT resolution is a lot coarser than used in the isolated pulsar case.  Therefore equal window sizes  (in bins) actually perform differently. 
  
• We have chosen a window size of 60 for historical reasons, this corresponds to  a size of 1 Hz.  In practical terms this sets a scale such that features broader than this are smoothed away but finer features remain.
• Note that the isolated search use a window size of 50, equivelent to 0.028 Hz.

At this point some output files are opened if defined at the top of the code by FILE_FMAX and FILE_FSTATS.  Respectively, these files output the maximum* Fstatistic + corresponding search parameters for each orbital search template within the search band and the resulting F statistics + search parameters for every single orbital + frequency filter searched.   Note that the * indicates that by maximum we mean maximum remaining F statistic once clusters of F statistics have been identified, recorded, and removed. 

Following this stage if the binary flag is not set then the code proceeds to setting up a grid of sky positions relavent to the isolated pulsar search.  If the binary flag is set then we skip this and move on to the following.

• ReadBinaryTemplateBank()
• Description : This function reads in a pre-generated file of orbital search templates and does some simple validation.
  
• Key Points
• Currently a single binary search template consists of the following quantities -
  
• The projected semi-major axis
• The orbital period
• The observed time of orbital periapse passage as measured in the SSB frame
• The orbital eccentricity
• The argument of periapse
• The header information consists of the parameters that were used to generate the template bank ie. parameter boundaries, observation time, mismatch, number of filters, maximum search frequency etc...

Here we enter the loop that cycles over search templates.  In the isolated pulsar case these search templates are sky positions, and in the binary search these are the orbital filters read into memory by the function ReadBinaryTemplateBank().

• CreateBinaryDemodParams()
• Description : This function is the equivelent of  CreateDemodParams() and calculates the input parameters neccessary for ComputeSkyBinary().  It also calls ComputeSkyBinary() to generate the input variables neccessary for LALDemod().
• Key Points : The difference between this and CreateDemodParams() is that a different data structure is filled in which has extra elements for containing the orbital parameters.

There is now a loop over spin down templates within which the value of the first order spin down parameter is simply incremented by a user defined input value.   For simplicity at present the value of the first order spin down parameter is set at its default value of 0.0 and is not incremented.  If the binary flag is set then we do not cycle over this loop, we pass through it only once with all spin down values set at zero.

• LALDemod()
• Description : This is the heart of the search code and is unchanged for the binary search.  It basically does the number crunching in the searches and will dominate all CPU time.

•  EstimateFLines()
• Description : This function takes the output from a single orbital/Sky position template and using the running median identifies "clusters" of F-statistic's above a user defined threshold.
• Key Points
• The window size we use here is hardcoded in (just as in the isolated search) at 100.
  
• The estimated maximum width of a true GW signal from a binary source for our chosen observation time has yet to be resolved.

If the user has specified an output file then at this point the results of the search over the current orbital filter are saved to file.  This includes the value of the detection statistic (which is in fact 2F X normalisation factor due to the descrete calculation of the running median) and the corresponding search parameters (orbital + frequency).

•  writeFLines()
• Description : This function takes the output of  EstimateFLines() and outputs the identified clusters to file.  It also removes the clusters from the data structures they were stored in leaving only non-cluster data remaining.

• writeFaFb()

• Description : This function writes to file the complex constituents Fa and Fb of the F statistics identified as parts of clusters.  This is done for later use in the Fstatistic shape test.

• EstimateSignalParameters()

• Description : This function takes the components Fa and Fb of the F statistics identified as parts of clusters and calculates estimated values of the parameters (polaristaion angle, inclination angle, and initial GW phase).  These are recorded for later use in the Fstatistic shape test.

• FreeMem()

• Description : This simply frees up the memory previously allocated in the code.