Loading [MathJax]/extensions/TeX/AMSsymbols.js
LALInspiral 5.0.3.1-ea7c608
All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Modules Pages
LALInspiralWave.c
Go to the documentation of this file.
1/*
2* Copyright (C) 2007 David Churches, Duncan Brown, Jolien Creighton, David McKechan, B.S. Sathyaprakash, Thomas Cokelaer, Riccardo Sturani, Laszlo Vereb
3*
4* This program is free software; you can redistribute it and/or modify
5* it under the terms of the GNU General Public License as published by
6* the Free Software Foundation; either version 2 of the License, or
7* (at your option) any later version.
8*
9* This program is distributed in the hope that it will be useful,
10* but WITHOUT ANY WARRANTY; without even the implied warranty of
11* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12* GNU General Public License for more details.
13*
14* You should have received a copy of the GNU General Public License
15* along with with program; see the file COPYING. If not, write to the
16* Free Software Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston,
17* MA 02110-1301 USA
18*/
19
20/**
21 * \author Churches, D. K and Sathyaprakash, B.S. and Cokelaer. T
22 * \file
23 * \ingroup LALInspiral_h
24 *
25 * \brief Interface routine needed to generate all waveforms in the \ref lalinspiral_inspiral package.
26 *
27 * To generate a waveform
28 * a user is noramlly required to (a) choose the binary parameters, starting frequency, number of
29 * bins of leading and trailing zero-padding, etc., in the structure <tt>InspiralTemplate params</tt>
30 * and (b) call the following three functions in the order given:
31 * LALInspiralParameterCalc(), LALInspiralWaveLength() and LALInspiralWave().
32 * Either a time- or a frequency-domain signal is returned depending upon the
33 * \c approximant requested (see Notes below).
34 *
35 * ### Prototypes ###
36 *
37 * <tt>LALInspiralWave()</tt>
38 * <ul>
39 * <li> \c signalvec: Output containing the inspiral waveform.</li>
40 * <li> \c params: Input containing binary chirp parameters.</li>
41 * </ul>
42 *
43 * <tt>LALInspiralWaveTemplates()</tt>
44 * <ul>
45 * <li> \c signalvec1: Output containing the 0-phase inspiral waveform.</li>
46 * <li> \c signalvec2: Output containing the \f$\pi/2\f$-phase inspiral waveform.</li>
47 * <li> \c params: Input containing binary chirp parameters.</li>
48 * </ul>
49 *
50 * ### Description ###
51 *
52 * The code LALInspiralWave() is the user interface to the inspiral codes. It takes from the user all
53 * the physical parameters which specify the binary, and calls the relevent wave generation function.
54 * Currently ten different approximants are fully implemented. These are #TaylorT1, #TaylorT2,
55 * #TaylorT3, #TaylorF1, #TaylorF2, #TaylorF2RedSpin, #PadeT1, #EOB, #BCV, #SpinTaylorT3, #PhenSpinTaylorRD.
56 * \c Taylor approximants can all be generated at seven different post-Newtonian orders,
57 * from Newtonian to 3.5 PN order, #PadeT1 exists at order 1.5PN and higher,
58 * #EOB at orders 2 and higher. #SpinTaylorT3 is implemented only at 2PN order
59 * by solving the evolution equations for the spin and orbital angular momenta and a
60 * time-domain phasing formula. Finally, PN order is undefined for #BCV.
61 * The approximant and the order are set up by the enums \c Approximant and \c LALPNOrder,
62 * respectively.
63 *
64 * The waveforms are all terminated one bin before the last stable orbit is reached.
65 * The last stable orbit corresponding to a given \c Approximant and \c LALPNOrder is
66 * defined as follows: For all \c Taylor approximants at orders 0PN, 1PN and 1.5PN
67 * \f$v_\textrm{lso}^2=1/6,\f$ and at 2PN, 2.5PN, 3PN and 3.5PN
68 * \f$v_\textrm{lso}^2 = x^\textrm{lso}_{T_4},\f$ where \f$x^\textrm{lso}_{T_4}\f$ is
69 * defined in \ref table_energy "this table". In the case of \c Pade approximant
70 * at 1.5PN order \f$v_\textrm{lso}^2=1/6,\f$ and at orders 2PN, 2.5PN, 3PN and 3.5PN
71 * \f$v_\textrm{lso}^2 = x^\textrm{lso}_{P_4},\f$ where \f$x^\textrm{lso}_{P_4}\f$ is
72 * defined in \ref table_energy "this table". In the case of #EOB approximant,
73 * defined only at orders greater than 2PN, the plunge waveform is
74 * terminated at the light-ring orbit defined by \eqref{eq_LightRing}.
75 *
76 * In the case of LALInspiralWaveTemplates() <tt>*signalvec1</tt>
77 * contains the `0-phase' inspiral template and <tt>*signalvec2</tt> contains
78 * a signal that is \f$\pi/2\f$ out of phase with respect to <tt>*signalvec1.</tt>
79 * Currently, a template pair is generated only for the following \c approximants:
80 * #TaylorT1, #TaylorT2, #TaylorT3, #PadeT1, #EOB.
81 *
82 * See the test codes for examples of how to generate different approximations.
83 *
84 * ### Algorithm ###
85 *
86 * Simple use of \c switch statement to access different PN approximations.
87 *
88 * ### Uses ###
89 *
90 * Depending on the user inputs one of the following functions is called:
91 * \code
92 * LALInspiralWave1()
93 * XLALInspiralWave2()
94 * LALInspiralWave3()
95 * XLALInspiralStationaryPhaseApprox1()
96 * XLALInspiralStationaryPhaseApprox2()
97 * LALEOBWaveform()
98 * LALBCVWaveform()
99 * LALInspiralSpinModulatedWave()
100 * \endcode
101 *
102 * ### Notes ###
103 *
104 * <ul>
105 * <li> A time-domain waveform is returned when the \c Approximant is one of
106 * #TaylorT1, #TaylorT2, #TaylorT3, #PadeT1, #EOB, #SpinTaylorT3, #PhenSpinTaylorRD, #SpinQuadTaylor
107 * </li><li> A frequency-domain waveform is returned when the \c Approximant is one of
108 * #TaylorF1, #TaylorF2, #TaylorF2RedSpin, #BCV.
109 * In these cases the code returns the real and imagninary parts of the
110 * Fourier domain signal in the convention of fftw. For a signal vector
111 * of length <tt>n=signalvec->length</tt> (\c n even):</li>
112 * <ul>
113 * <li> <tt>signalvec->data[0]</tt> is the \e real 0th frequency component of the Fourier transform.</li>
114 * <li> <tt>signalvec->data[n/2]</tt> is the \e real Nyquist frequency component of the Fourier transform.</li>
115 * <li> <tt>signalvec->data[k]</tt> and <tt>signalvec->data[n-k],</tt> for <tt>k=1,..., n/2-1,</tt> are
116 * the real and imaginary parts of the Fourier transform at a frequency \f$k\Delta f=k/T,\f$ \f$T\f$ being
117 * the duration of the signal and \f$\Delta f=1/T\f$ is the frequency resolution.</li>
118 * </ul>
119 * </ul>
120 *
121 */
122
123#include <ctype.h>
124#include <lal/LALInspiral.h>
125#include <lal/LALNoiseModels.h>
126#include <lal/LALStdlib.h>
127#include <lal/GeneratePPNInspiral.h>
128#include <lal/LALSQTPNWaveformInterface.h>
129#include <lal/TimeSeries.h>
130#include <lal/LALString.h>
131#include <lal/LALDict.h>
132
133/**
134 * Generate the plus and cross polarizations for a waveform
135 * form a row of the sim_inspiral table.
136 *
137 * Parses a row from the sim_inspiral table and passes the appropriate members
138 * to XLALSimInspiralChooseWaveform().
139 *
140 * FIXME: this should eventually be moved to lalsimulation
141 * along with the appropriate string parsing functions
142 */
144 REAL8TimeSeries **hplus, /**< +-polarization waveform */
145 REAL8TimeSeries **hcross, /**< x-polarization waveform */
146 SimInspiralTable *thisRow, /**< row from the sim_inspiral table containing waveform parameters */
147 REAL8 deltaT /**< time step */
148 )
149{
150 int ret;
151 LALPNOrder order;
154
155 REAL8 phi0 = thisRow->coa_phase;
156 REAL8 m1 = thisRow->mass1 * LAL_MSUN_SI;
157 REAL8 m2 = thisRow->mass2 * LAL_MSUN_SI;
158 REAL8 S1x = thisRow->spin1x;
159 REAL8 S1y = thisRow->spin1y;
160 REAL8 S1z = thisRow->spin1z;
161 REAL8 S2x = thisRow->spin2x;
162 REAL8 S2y = thisRow->spin2y;
163 REAL8 S2z = thisRow->spin2z;
164 REAL8 f_min = thisRow->f_lower;
165 REAL8 f_ref = 0.;
166 REAL8 r = thisRow->distance * LAL_PC_SI * 1e6;
167 REAL8 i = thisRow->inclination;
168 LALDict *LALpars=XLALCreateDict();
170
171 /* get approximant */
173 if ( (int) approximant == XLAL_FAILURE)
175
176 /* get phase PN order; this is an enum such that the value is twice the PN order */
178 if ( (int) order == XLAL_FAILURE)
181
182 /* get taper option */
184 if ( (int) taper == XLAL_FAILURE)
186
187 /* generate +,x waveforms */
188 /* special case for NR waveforms */
189 switch(approximant)
190 {
191 case NumRelNinja2:
192 if (XLALNRInjectionFromSimInspiral(hplus, hcross, thisRow, deltaT) == XLAL_FAILURE) {
193 XLALDestroyDict(LALpars);
195 }
196 XLALDestroyDict(LALpars);
197 break;
198
199 default:
200 ret = XLALSimInspiralChooseTDWaveform(hplus, hcross, m1, m2,
201 S1x, S1y, S1z, S2x, S2y, S2z, r, i,
202 phi0, 0., 0., 0., deltaT, f_min, f_ref,
203 LALpars, approximant);
204 XLALDestroyDict(LALpars);
205 if( ret == XLAL_FAILURE )
207 }
208
209 /* taper the waveforms */
210 if (XLALSimInspiralREAL8WaveTaper((*hplus)->data, taper) == XLAL_FAILURE)
212
213 if (XLALSimInspiralREAL8WaveTaper((*hcross)->data, taper) == XLAL_FAILURE)
215
216 return XLAL_SUCCESS;
217}
218
219/**
220 * Generate the plus and cross polarizations for a conditioned waveform
221 * form a row of the sim_inspiral table.
222 *
223 * Parses a row from the sim_inspiral table and passes the appropriate members
224 * to XLALSimInspiralTD().
225 *
226 */
228 REAL8TimeSeries **hplus, /**< +-polarization waveform */
229 REAL8TimeSeries **hcross, /**< x-polarization waveform */
230 SimInspiralTable *thisRow, /**< row from the sim_inspiral table containing waveform parameters */
231 REAL8 deltaT /**< time step (s) */
232 )
233{
234 int ret;
235 LALPNOrder order;
237 REAL8 phi0 = thisRow->coa_phase;
238 REAL8 m1 = thisRow->mass1 * LAL_MSUN_SI;
239 REAL8 m2 = thisRow->mass2 * LAL_MSUN_SI;
240 REAL8 S1x = thisRow->spin1x;
241 REAL8 S1y = thisRow->spin1y;
242 REAL8 S1z = thisRow->spin1z;
243 REAL8 S2x = thisRow->spin2x;
244 REAL8 S2y = thisRow->spin2y;
245 REAL8 S2z = thisRow->spin2z;
246 REAL8 f_min = thisRow->f_lower;
247 REAL8 f_ref = thisRow->f_final;
248 REAL8 distance = thisRow->distance * LAL_PC_SI * 1e6;
249 REAL8 inclination = thisRow->inclination;
250 LALDict *params = XLALCreateDict();
252
253 /* get approximant */
255 if ( (int) approximant == XLAL_FAILURE)
257
258 if (approximant == NR_hdf5) {
259 char *filepath = thisRow->numrel_data;
261 }
262
263 /* get phase PN order; this is an enum such that the value is twice the PN order */
265 if ( (int) order == XLAL_FAILURE)
268
269 /* note: the condition waveform already does tapering... ignore any request to do so get taper option */
270 /* taper = XLALGetTaperFromString(thisRow->taper); */
271
272 /* get extra parameters from numrel_data column */
273 /* format is: "param1=3.14,param2:INT4=7,param3:string=hello" */
274 if (*thisRow->numrel_data && approximant != NR_hdf5) {
275 char *str = XLALStringDuplicate(thisRow->numrel_data);
276 char *head = str;
277 char *tok;
278 int i = 0;
279 /* remove whitespace from string */
280 for (int j=0; str[j]; ++j)
281 if (!isspace(str[j]))
282 str[i++] = str[j];
283 str[i] = '\0';
284 while ((tok = XLALStringToken(&str, ",", 0))) {
285 char *key = tok;
286 char *val = strchr(tok, '=');
287 char *type;
288 if (val == NULL) {
290 XLALFree(head);
291 XLAL_ERROR(XLAL_ENAME, "could not parse extra params from numrel_data column `%s'", thisRow->numrel_data);
292 }
293 *val++ = 0;
294 type = strchr(key, ':');
295 if (type)
296 *type++ = 0;
297 if (type == NULL || XLALStringCaseCompare(type, "REAL8") == 0)
299 else if (XLALStringCaseCompare(type, "INT4") == 0)
301 else if (XLALStringCaseCompare(type, "STRING") == 0)
303 else {
305 XLALFree(head);
306 XLAL_ERROR(XLAL_ENAME, "could not parse extra params from numrel_data column `%s'", thisRow->numrel_data);
307 }
308 }
309 XLALFree(head);
310 }
311
312 /* generate +,x waveforms */
313 ret = XLALSimInspiralTD(hplus, hcross, m1, m2, S1x, S1y, S1z, S2x, S2y, S2z, distance, inclination, phi0, 0.0, 0.0, 0.0, deltaT, f_min, f_ref, params, approximant);
314
316 if( ret == XLAL_FAILURE )
318
319 return XLAL_SUCCESS;
320}
321
322/**
323 * Generate the plus and cross polarizations for a waveform
324 * form a row of the InspiralTemplate structure.
325 *
326 * Parses the InspiralTemplate stucture and passes the appropriate members
327 * to XLALSimInspiralChooseWaveform().
328 */
329int
331 REAL8TimeSeries **hplus, /**< +-polarization waveform */
332 REAL8TimeSeries **hcross, /**< x-polarization waveform */
333 InspiralTemplate *params /**< stucture containing waveform parameters */
334 )
335{
336 int ret;
337 REAL8 deltaT = 1./params->tSampling;
338 REAL8 phi0 = params->startPhase; /* startPhase is used as the peak amplitude phase here */
339 REAL8 m1 = params->mass1 * LAL_MSUN_SI;
340 REAL8 m2 = params->mass2 * LAL_MSUN_SI;
341 REAL8 S1x = params->spin1[0];
342 REAL8 S1y = params->spin1[1];
343 REAL8 S1z = params->spin1[2];
344 REAL8 S2x = params->spin2[0];
345 REAL8 S2y = params->spin2[1];
346 REAL8 S2z = params->spin2[2];
347 REAL8 f_min = params->fLower;
348 REAL8 f_ref = 0.;
349 /* Value of 'distance' fed to lalsim is conventional to obtain a correct template norm */
352 LALDict *LALpars=XLALCreateDict();
355 Approximant approximant = params->approximant;
356
357 /* generate +,x waveforms */
358 ret = XLALSimInspiralChooseTDWaveform(hplus, hcross, m1, m2,
359 S1x, S1y, S1z, S2x, S2y, S2z, r, i,
360 phi0, 0., 0., 0., deltaT, f_min, f_ref,
361 LALpars, approximant);
362 XLALDestroyDict(LALpars);
363
364 if( ret == XLAL_FAILURE)
366
367 return XLAL_SUCCESS;
368}
369
370
371void
374 REAL4Vector *signalvec,
376 )
377{
378
381
382 ASSERT (signalvec, status, LALINSPIRALH_ENULL, LALINSPIRALH_MSGENULL);
383 ASSERT (signalvec->length >= 2, status, LALINSPIRALH_ESIZE, LALINSPIRALH_MSGESIZE);
384 ASSERT (signalvec->data, status, LALINSPIRALH_ENULL, LALINSPIRALH_MSGENULL);
385 ASSERT (params, status, LALINSPIRALH_ENULL, LALINSPIRALH_MSGENULL);
386
387
388 ASSERT((INT4)params->approximant >= 0, status, LALINSPIRALH_ESIZE, LALINSPIRALH_MSGESIZE);
389 ASSERT((INT4)params->approximant < NumApproximants, status, LALINSPIRALH_ESIZE, LALINSPIRALH_MSGESIZE);
391 status, LALINSPIRALH_EORDER, LALINSPIRALH_MSGEORDER);
392
394 {
395 REAL8TimeSeries *hplus = NULL;
396 REAL8TimeSeries *hcross = NULL;
397 unsigned int idx;
398
399 /* generate hplus and hcross */
402
403 /* check length of waveform compared to signalvec */
404 if (hplus->data->length > signalvec->length)
405 {
408 XLALPrintError("XLALSimInspiralChooseWaveformFromInspiralTemplate generated a waveform longer than output vector.\n");
409 ABORT(status, LALINSPIRALH_EVECTOR, LALINSPIRALH_MSGEVECTOR);
410 }
411
412 /* initialize template vector to zero since hplus may be shorter */
413 memset(signalvec->data, 0.0, signalvec->length * sizeof(signalvec->data[0]));
414
415 /* convert REAL8 hplus waveform to REAL4 and store it in signalvec->data */
416 for(idx = 0; idx < signalvec->length && idx < hplus->data->length ; idx++)
417 {
418 signalvec->data[idx] = (REAL4) hplus->data->data[idx];
419 }
420
421 params->tC = - ( ((REAL8)hplus->epoch.gpsSeconds) + 1.e-9*((REAL8)hplus->epoch.gpsNanoSeconds) );
422
423 /* free hplus and hcross */
426 }
427 else switch (params->approximant)
428 {
429 case TaylorT1:
430 case PadeT1:
432 break;
433 case TaylorT2:
435 break;
436 case TaylorT3:
438 break;
439 case EOB:
440 case EOBNR:
442 break;
443 case EOBNRv2:
444 case EOBNRv2HM:
446 break;
447 case IMRPhenomA:
448 case IMRPhenomB:
450 break;
451 case IMRPhenomFA:
453 break;
454 case IMRPhenomFB:
456 break;
457 case BCV:
458 LALBCVWaveform(status->statusPtr, signalvec, params);
460 break;
461 case BCVSpin:
462 LALBCVSpinWaveform(status->statusPtr, signalvec, params);
464 break;
465 case TaylorF1:
467 break;
468 case TaylorF2:
469 case FindChirpSP:
471 break;
472 case TaylorF2RedSpin:
474 break;
475 case PadeF1:
476 ABORT(status, LALINSPIRALH_ECHOICE, LALINSPIRALH_MSGECHOICE);
477 break;
478 case SpinTaylorT3:
479 LALInspiralSpinModulatedWave(status->statusPtr, signalvec, params);
481 break;
482 case SpinTaylor:
483 /* this generate the h+ waveform whereas the one above takes into
484 * account h+, hx and orientation of the system*/
485 LALSTPNWaveform(status->statusPtr, signalvec, params);
487 break;
490 break;
491 case PhenSpinTaylorRD:
493 break;
494 case SpinQuadTaylor:
497 break;
498 case AmpCorPPN:
499 LALInspiralAmplitudeCorrectedWave(status->statusPtr, signalvec, params);
501 break;
502 case Eccentricity:
503 LALInspiralEccentricity(status->statusPtr, signalvec, params);
505 break;
506 case TaylorEt:
508 break;
509 case TaylorT4:
510 LALTaylorT4Waveform(status->statusPtr, signalvec, params);
512 break;
513 case TaylorN:
515 break;
516 default:
517 ABORT( status, LALINSPIRALH_ESWITCH, LALINSPIRALH_MSGESWITCH );
518 }
519
521 RETURN (status);
522}
523
524void
527 REAL4Vector *signalvec1,
528 REAL4Vector *signalvec2,
530 )
531{
532
535
536 ASSERT (signalvec1, status, LALINSPIRALH_ENULL, LALINSPIRALH_MSGENULL);
537 ASSERT (signalvec2, status, LALINSPIRALH_ENULL, LALINSPIRALH_MSGENULL);
538 ASSERT (signalvec1->length >= 2, status, LALINSPIRALH_ESIZE, LALINSPIRALH_MSGESIZE);
539 ASSERT (signalvec2->length >= 2, status, LALINSPIRALH_ESIZE, LALINSPIRALH_MSGESIZE);
540 ASSERT (signalvec1->data, status, LALINSPIRALH_ENULL, LALINSPIRALH_MSGENULL);
541 ASSERT (signalvec2->data, status, LALINSPIRALH_ENULL, LALINSPIRALH_MSGENULL);
542 ASSERT (params, status, LALINSPIRALH_ENULL, LALINSPIRALH_MSGENULL);
543
544
545 ASSERT((INT4)params->approximant >= 0, status, LALINSPIRALH_ESIZE, LALINSPIRALH_MSGESIZE);
546 ASSERT((INT4)params->approximant < NumApproximants, status, LALINSPIRALH_ESIZE, LALINSPIRALH_MSGESIZE);
548 status, LALINSPIRALH_EORDER, LALINSPIRALH_MSGEORDER);
549
551 {
552 REAL8TimeSeries *hplus = NULL;
553 REAL8TimeSeries *hcross = NULL;
554 unsigned int idx;
555
556 /* generate hplus and hcross */
559
560 /* check length of waveform compared to signalvec */
561 if (hplus->data->length > signalvec1->length)
562 {
565 XLALPrintError("XLALSimInspiralChooseWaveformFromInspiralTemplate generated a waveform longer than output vector.\n");
566 ABORT(status, LALINSPIRALH_EVECTOR, LALINSPIRALH_MSGEVECTOR);
567 }
568
569 /* initialize template vectors to zero since hplus and hcross may be shorter */
570 memset(signalvec1->data, 0.0, signalvec1->length * sizeof(signalvec1->data[0]));
571 memset(signalvec2->data, 0.0, signalvec2->length * sizeof(signalvec2->data[0]));
572
573 /* convert REAL8 hplus and hcross waveforms to REAL4 and store them in signalvec[1,2]->data */
574 for(idx = 0; idx < signalvec1->length; idx++)
575 {
576 signalvec1->data[idx] = (REAL4) hplus->data->data[idx];
577 signalvec2->data[idx] = (REAL4) hcross->data->data[idx];
578 }
579
580 params->tC = - ( ((REAL8)hplus->epoch.gpsSeconds) + 1.e-9*((REAL8)hplus->epoch.gpsNanoSeconds) );
581
582 /* free hplus and hcross */
585 }
586 else switch (params->approximant)
587 {
588 case TaylorT1:
589 case PadeT1:
590 if (XLALInspiralWave1Templates(signalvec1, signalvec2, params) == XLAL_FAILURE) ABORTXLAL(status);
591 break;
592 case TaylorT2:
593 if (XLALInspiralWave2Templates(signalvec1, signalvec2, params) == XLAL_FAILURE) ABORTXLAL(status);
594 break;
595 case TaylorT3:
596 if (XLALInspiralWave3Templates(signalvec1, signalvec2, params) == XLAL_FAILURE) ABORTXLAL(status);
597 break;
598 case TaylorEt:
599 if (XLALTaylorEtWaveformTemplates(signalvec1, signalvec2, params) == XLAL_FAILURE) ABORTXLAL(status);
600 break;
601 case EOB:
602 case EOBNR:
603 if (XLALEOBWaveformTemplates(signalvec1, signalvec2, params) == XLAL_FAILURE) ABORTXLAL(status);
604 break;
605 case EOBNRv2:
606 case EOBNRv2HM:
607 if (XLALEOBPPWaveformTemplates(signalvec1, signalvec2, params) == XLAL_FAILURE) ABORTXLAL(status);
608 break;
609 case IMRPhenomA:
610 case IMRPhenomB:
611 if (XLALBBHPhenWaveTimeDomTemplates(signalvec1, signalvec2, params) == XLAL_FAILURE) ABORTXLAL(status);
612 break;
613 case IMRPhenomFA:
614 if (XLALBBHPhenWaveAFreqDomTemplates(signalvec1, signalvec2, params) == XLAL_FAILURE) ABORTXLAL(status);
615 break;
616 case IMRPhenomFB:
617 if (XLALBBHPhenWaveBFreqDomTemplates(signalvec1, signalvec2, params) == XLAL_FAILURE) ABORTXLAL(status);
618 break;
619 case TaylorF2RedSpin:
620 if (XLALTaylorF2ReducedSpinTemplates(signalvec1, signalvec2, params) == XLAL_FAILURE) ABORTXLAL(status);
621 break;
622 case TaylorF1:
623 case TaylorF2:
624 case FindChirpSP:
625 case PadeF1:
626 case BCV:
627 case BCVSpin:
628 ABORT(status, LALINSPIRALH_ECHOICE, LALINSPIRALH_MSGECHOICE);
629 break;
630 case SpinTaylor:
631 LALSTPNWaveformTemplates(status->statusPtr, signalvec1, signalvec2, params);
633 break;
636 break;
637 case PhenSpinTaylorRD:
638 if (XLALPSpinInspiralRDTemplates(signalvec1, signalvec2, params) == XLAL_FAILURE) ABORTXLAL(status);
639 break;
640 case SpinQuadTaylor:
641 LALSTPNWaveformTemplates(status->statusPtr, signalvec1, signalvec2, params);
642 break;
643 case AmpCorPPN:
644 LALInspiralAmplitudeCorrectedWaveTemplates(status->statusPtr, signalvec1, signalvec2, params);
646 break;
647 case Eccentricity:
648 LALInspiralEccentricityTemplates(status->statusPtr, signalvec1, signalvec2, params);
650 break;
651 default:
652 ABORT( status, LALINSPIRALH_ESWITCH, LALINSPIRALH_MSGESWITCH );
653
654 }
655
657 RETURN (status);
658}
659
660void
663 CoherentGW *waveform,
664 InspiralTemplate *inspiralParams,
665 PPNParamStruc *ppnParams)
666{
667
670
671 ASSERT (inspiralParams, status, LALINSPIRALH_ENULL, LALINSPIRALH_MSGENULL);
672 ASSERT (ppnParams, status, LALINSPIRALH_ENULL, LALINSPIRALH_MSGENULL);
673
674
676 {
677 REAL8TimeSeries *hplus = NULL;
678 REAL8TimeSeries *hcross = NULL;
679 unsigned int idx;
680
681 /* generate hplus and hcross */
682 if (XLALSimInspiralChooseWaveformFromInspiralTemplate(&hplus, &hcross, inspiralParams) == XLAL_FAILURE)
684
685 /* allocate the waveform data stucture h in the CoherentGW structure */
686 waveform->h = (REAL4TimeVectorSeries *) LALMalloc( sizeof(REAL4TimeVectorSeries) );
687 if ( waveform->h == NULL )
688 {
691 ABORT(status, LALINSPIRALH_EMEM, LALINSPIRALH_MSGEMEM);
692 }
693 memset( waveform->h, 0, sizeof(REAL4TimeVectorSeries) );
694
695 /* populate the waveform data stucture h in the CoherentGW structure */
696 waveform->h->data = XLALCreateREAL4VectorSequence(hplus->data->length, 2);
697 if (waveform->h == NULL)
698 {
699 CHAR warnMsg[1024];
702 LALFree(waveform->h);
703 snprintf( warnMsg, XLAL_NUM_ELEM(warnMsg),
704 "Memory allocation error when allocating CoherentGW REAL4VectorSequence.\n");
705 LALInfo( status, warnMsg );
706 ABORT( status, LALINSPIRALH_EMEM, LALINSPIRALH_MSGEMEM );
707 }
708 waveform->h->f0 = hplus->f0;
709 waveform->h->deltaT = hplus->deltaT;
710 waveform->h->epoch = hplus->epoch;
711 waveform->h->sampleUnits = hplus->sampleUnits;
712 waveform->position = ppnParams->position;
713 waveform->psi = ppnParams->psi;
714
715 /* convert REAL8 hplus and hcross waveforms to REAL4 and store them in waveform->h->data */
716 for(idx = 0; idx < waveform->h->data->length; idx++)
717 {
718 waveform->h->data->data[idx * 2] = (REAL4) hplus->data->data[idx];
719 waveform->h->data->data[idx * 2 + 1] = (REAL4) hcross->data->data[idx];
720 }
721
722 /* free hplus and hcross */
725 }
726 else switch (inspiralParams->approximant)
727 {
728 case TaylorT1:
729 case PadeT1:
730 if (XLALInspiralWave1ForInjection(waveform, inspiralParams, ppnParams) == XLAL_FAILURE) ABORTXLAL(status);
731 break;
732 case TaylorT2:
733 if (XLALInspiralWave2ForInjection(waveform, inspiralParams, ppnParams) == XLAL_FAILURE) ABORTXLAL(status);
734 break;
735 case TaylorT3:
736 if (XLALInspiralWave3ForInjection(waveform, inspiralParams, ppnParams) == XLAL_FAILURE) ABORTXLAL(status);
737 break;
738 case EOB:
739 case EOBNR:
740 if (XLALEOBWaveformForInjection(waveform, inspiralParams, ppnParams) == XLAL_FAILURE) ABORTXLAL(status);
741 break;
742 case EOBNRv2:
743 case EOBNRv2HM:
744 if (XLALEOBPPWaveformForInjection(waveform, inspiralParams, ppnParams) == XLAL_FAILURE) ABORTXLAL(status);
745 break;
746 case IMRPhenomA:
747 case IMRPhenomB:
748 if (XLALBBHPhenWaveTimeDomForInjection(waveform, inspiralParams, ppnParams) == XLAL_FAILURE) ABORTXLAL(status);
749 break;
750 case BCV:
751 /*LALBCVWaveformForInjection(status->statusPtr, waveform, inspiralParams, ppnParams);*/
752 case BCVSpin:
753 case TaylorF1:
754 case TaylorF2:
755 case TaylorF2RedSpin:
756 case FindChirpSP:
757 case PadeF1:
758 ABORT(status, LALINSPIRALH_ECHOICE, LALINSPIRALH_MSGECHOICE);
759 break;
761 if (XLALSTPNFramelessWaveformForInjection(waveform, inspiralParams, ppnParams) == XLAL_FAILURE) ABORTXLAL(status);
762 break;
763 case SpinTaylor:
764 LALSTPNWaveformForInjection(status->statusPtr, waveform, inspiralParams, ppnParams);
766 break;
767 case PhenSpinTaylorRD:
768 if (XLALPSpinInspiralRDForInjection(waveform, inspiralParams, ppnParams) == XLAL_FAILURE) ABORTXLAL(status);
769 break;
770 case SpinQuadTaylor:
771 if (XLALSQTPNWaveformForInjection(waveform, inspiralParams, ppnParams) == XLAL_FAILURE) ABORTXLAL(status);
772 break;
773 case AmpCorPPN:
774 LALInspiralAmplitudeCorrectedWaveForInjection(status->statusPtr, waveform, inspiralParams, ppnParams);
776 break;
777 case TaylorT4:
778 LALTaylorT4WaveformForInjection(status->statusPtr, waveform, inspiralParams, ppnParams);
780 break;
781 default:
782 ABORT( status, LALINSPIRALH_ESWITCH, LALINSPIRALH_MSGESWITCH );
783
784 }
785
787 RETURN (status);
788}
void LALBCVWaveform(LALStatus *status, REAL4Vector *signalvec, InspiralTemplate *params)
void LALBCVSpinWaveform(LALStatus *status, REAL4Vector *signalvec, InspiralTemplate *params)
void XLALDestroyDict(LALDict *dict)
LALDict * XLALCreateDict(void)
int XLALDictInsertStringValue(LALDict *dict, const char *key, const char *string)
int XLALDictInsertINT4Value(LALDict *dict, const char *key, INT4 value)
int XLALDictInsertREAL8Value(LALDict *dict, const char *key, REAL8 value)
int XLALEOBPPWaveformForInjection(CoherentGW *waveform, InspiralTemplate *params, PPNParamStruc *ppnParams)
int XLALEOBPPWaveformTemplates(REAL4Vector *signalvec1, REAL4Vector *signalvec2, InspiralTemplate *params)
int XLALEOBPPWaveform(REAL4Vector *signalvec, InspiralTemplate *params)
int XLALEOBWaveform(REAL4Vector *signalvec, InspiralTemplate *params)
int XLALEOBWaveformForInjection(CoherentGW *waveform, InspiralTemplate *params, PPNParamStruc *ppnParams)
int XLALEOBWaveformTemplates(REAL4Vector *signalvec1, REAL4Vector *signalvec2, InspiralTemplate *params)
void LALSTPNWaveform(LALStatus *status, REAL4Vector *signalvec, InspiralTemplate *params)
int XLALInspiralWave2ForInjection(CoherentGW *waveform, InspiralTemplate *params, PPNParamStruc *ppnParams)
int XLALBBHPhenWaveTimeDomTemplates(REAL4Vector *signalvec1, REAL4Vector *signalvec2, InspiralTemplate *insp_template)
int XLALInspiralStationaryPhaseApprox2(REAL4Vector *signalvec, InspiralTemplate *params)
int XLALSTPNFramelessWaveform(REAL4Vector *signalvec, InspiralTemplate *params)
void LALInspiralEccentricity(LALStatus *status, REAL4Vector *signalvec, InspiralTemplate *params)
void LALInspiralAmplitudeCorrectedWaveTemplates(LALStatus *status, REAL4Vector *filter1, REAL4Vector *filter2, InspiralTemplate *params)
int XLALInspiralStationaryPhaseApprox1(REAL4Vector *signalvec, InspiralTemplate *params)
void LALInspiralSpinModulatedWave(LALStatus *status, REAL4Vector *signalvec, InspiralTemplate *in)
int XLALInspiralWave1(REAL4Vector *signalvec, InspiralTemplate *params)
int XLALInspiralWave1Templates(REAL4Vector *signalvec1, REAL4Vector *signalvec2, InspiralTemplate *params)
void LALInspiralAmplitudeCorrectedWaveForInjection(LALStatus *status, CoherentGW *waveform, InspiralTemplate *params, PPNParamStruc *ppnParams)
int XLALBBHPhenWaveAFreqDomTemplates(REAL4Vector *signalvec1, REAL4Vector *signalvec2, InspiralTemplate *params)
int XLALInspiralWave2Templates(REAL4Vector *signalvec1, REAL4Vector *signalvec2, InspiralTemplate *params)
int XLALTaylorNWaveform(REAL4Vector *signalvec, InspiralTemplate *params)
int XLALInspiralWave1ForInjection(CoherentGW *waveform, InspiralTemplate *params, PPNParamStruc *ppnParams)
int XLALBBHPhenWaveBFreqDom(REAL4Vector *signalvec, InspiralTemplate *params)
int XLALInspiralWave2(REAL4Vector *signalvec, InspiralTemplate *params)
void LALInspiralAmplitudeCorrectedWave(LALStatus *status, REAL4Vector *signalvec, InspiralTemplate *params)
int XLALNRInjectionFromSimInspiral(REAL8TimeSeries **hplus, REAL8TimeSeries **hcross, SimInspiralTable *thisRow, REAL8 deltaT)
int XLALBBHPhenWaveBFreqDomTemplates(REAL4Vector *signalvec1, REAL4Vector *signalvec2, InspiralTemplate *params)
int XLALSTPNFramelessWaveformForInjection(CoherentGW *waveform, InspiralTemplate *params, PPNParamStruc *ppnParams)
void LALTaylorT4Waveform(LALStatus *status, REAL4Vector *signalvec, InspiralTemplate *params)
int XLALInspiralWave3Templates(REAL4Vector *signalvec1, REAL4Vector *signalvec2, InspiralTemplate *params)
int XLALTaylorF2ReducedSpin(REAL4Vector *signalvec, InspiralTemplate *params)
Generate the "reduced-spin templates" proposed in http://arxiv.org/abs/1107.1267.
int XLALBBHPhenWaveTimeDomForInjection(CoherentGW *waveform, InspiralTemplate *params, PPNParamStruc *ppnParams)
int XLALTaylorF2ReducedSpinTemplates(REAL4Vector *signalvec1, REAL4Vector *signalvec2, InspiralTemplate *params)
Generate two orthogonal "reduced-spin" templates.
int XLALBBHPhenWaveAFreqDom(REAL4Vector *signalvec, InspiralTemplate *params)
int XLALInspiralWave3(REAL4Vector *signalvec, InspiralTemplate *params)
void LALSTPNWaveformTemplates(LALStatus *status, REAL4Vector *signalvec1, REAL4Vector *signalvec2, InspiralTemplate *params)
int XLALTaylorEtWaveform(REAL4Vector *signalvec, InspiralTemplate *params)
int XLALInspiralWave3ForInjection(CoherentGW *waveform, InspiralTemplate *params, PPNParamStruc *ppnParams)
int XLALBBHPhenWaveTimeDom(REAL4Vector *signalvec, InspiralTemplate *insp_template)
void LALSTPNWaveformForInjection(LALStatus *status, CoherentGW *waveform, InspiralTemplate *params, PPNParamStruc *ppnParams)
void LALInspiralEccentricityTemplates(LALStatus *status, REAL4Vector *signalvec1, REAL4Vector *signalvec2, InspiralTemplate *params)
int XLALSTPNFramelessWaveformTemplates(REAL4Vector *signalvec1, REAL4Vector *signalvec2, InspiralTemplate *params)
void LALTaylorT4WaveformForInjection(LALStatus *status, CoherentGW *waveform, InspiralTemplate *params, PPNParamStruc *ppnParams)
int XLALTaylorEtWaveformTemplates(REAL4Vector *signalvec1, REAL4Vector *signalvec2, InspiralTemplate *params)
int XLALInspiralTDWaveformFromSimInspiral(REAL8TimeSeries **hplus, REAL8TimeSeries **hcross, SimInspiralTable *thisRow, REAL8 deltaT)
Generate the plus and cross polarizations for a conditioned waveform form a row of the sim_inspiral t...
int XLALSimInspiralChooseWaveformFromInspiralTemplate(REAL8TimeSeries **hplus, REAL8TimeSeries **hcross, InspiralTemplate *params)
Generate the plus and cross polarizations for a waveform form a row of the InspiralTemplate structure...
void LALInspiralWaveTemplates(LALStatus *status, REAL4Vector *signalvec1, REAL4Vector *signalvec2, InspiralTemplate *params)
void LALInspiralWaveForInjection(LALStatus *status, CoherentGW *waveform, InspiralTemplate *inspiralParams, PPNParamStruc *ppnParams)
int XLALSimInspiralChooseWaveformFromSimInspiral(REAL8TimeSeries **hplus, REAL8TimeSeries **hcross, SimInspiralTable *thisRow, REAL8 deltaT)
Generate the plus and cross polarizations for a waveform form a row of the sim_inspiral table.
void LALInspiralWave(LALStatus *status, REAL4Vector *signalvec, InspiralTemplate *params)
#define LALMalloc(n)
#define LALFree(p)
int XLALSimInspiralGetTaperFromString(const char *string)
int XLALSimInspiralGetPNOrderFromString(const char *waveform)
int XLALSimInspiralImplementedTDApproximants(Approximant approximant)
int XLALSimInspiralGetApproximantFromString(const char *waveform)
int XLALSimInspiralWaveformParamsInsertPNPhaseOrder(LALDict *params, INT4 value)
int XLALSimInspiralWaveformParamsInsertNumRelData(LALDict *params, const char *value)
int XLALSimInspiralWaveformParamsInsertPNAmplitudeOrder(LALDict *params, INT4 value)
#define ABORT(statusptr, code, mesg)
#define CHECKSTATUSPTR(statusptr)
#define ATTATCHSTATUSPTR(statusptr)
#define ASSERT(assertion, statusptr, code, mesg)
#define DETATCHSTATUSPTR(statusptr)
#define INITSTATUS(statusptr)
#define RETURN(statusptr)
#define ABORTXLAL(sp)
double i
REAL4VectorSequence * XLALCreateREAL4VectorSequence(UINT4 length, UINT4 veclen)
#define LAL_MSUN_SI
#define LAL_PC_SI
#define XLAL_NUM_ELEM(x)
double REAL8
char CHAR
uint32_t UINT4
int32_t INT4
float REAL4
int LALInfo(LALStatus *status, const char *info)
#define LALINSPIRALH_ECHOICE
Invalid choice for an input parameter.
Definition: LALInspiral.h:60
#define LALINSPIRALH_EVECTOR
Attempting to write beyond the end of vector.
Definition: LALInspiral.h:71
#define LALINSPIRALH_EMEM
Memory allocation error.
Definition: LALInspiral.h:57
#define LALINSPIRALH_ENULL
Arguments contained an unexpected null pointer.
Definition: LALInspiral.h:56
#define LALINSPIRALH_EORDER
unknown order specified
Definition: LALInspiral.h:61
#define LALINSPIRALH_ESWITCH
Unknown case in switch.
Definition: LALInspiral.h:75
#define LALINSPIRALH_ESIZE
Invalid input range.
Definition: LALInspiral.h:59
void XLALFree(void *p)
int XLALPSpinInspiralRDTemplates(REAL4Vector *signalvec1, REAL4Vector *signalvec2, InspiralTemplate *params)
Function to produce waveform templates.
int XLALPSpinInspiralRD(REAL4Vector *signalvec, InspiralTemplate *params)
int XLALPSpinInspiralRDForInjection(CoherentGW *waveform, InspiralTemplate *params, PPNParamStruc *ppnParams)
Function Module to produce injection waveforms.
int XLALSQTPNWaveform(REAL4Vector *signalvec, InspiralTemplate *params)
The function returns the generated waveform.
int XLALSQTPNWaveformForInjection(CoherentGW *waveform, InspiralTemplate *params, PPNParamStruc *ppnParams)
The function returns the generated waveform for injection.
int XLALSimInspiralChooseTDWaveform(REAL8TimeSeries **hplus, REAL8TimeSeries **hcross, const REAL8 m1, const REAL8 m2, const REAL8 S1x, const REAL8 S1y, const REAL8 S1z, const REAL8 S2x, const REAL8 S2y, const REAL8 S2z, const REAL8 distance, const REAL8 inclination, const REAL8 phiRef, const REAL8 longAscNodes, const REAL8 eccentricity, const REAL8 UNUSED meanPerAno, const REAL8 deltaT, const REAL8 f_min, REAL8 f_ref, LALDict *params, const Approximant approximant)
int XLALSimInspiralTD(REAL8TimeSeries **hplus, REAL8TimeSeries **hcross, REAL8 m1, REAL8 m2, REAL8 S1x, REAL8 S1y, REAL8 S1z, REAL8 S2x, REAL8 S2y, REAL8 S2z, REAL8 distance, REAL8 inclination, REAL8 phiRef, REAL8 longAscNodes, REAL8 eccentricity, REAL8 meanPerAno, REAL8 deltaT, REAL8 f_min, REAL8 f_ref, LALDict *LALparams, Approximant approximant)
LALSimInspiralApplyTaper
Approximant
LALPNOrder
EOB
PadeT1
IMRPhenomFB
NumApproximants
EOBNRv2HM
TaylorEt
BCVSpin
EOBNRv2
TaylorN
IMRPhenomFA
SpinQuadTaylor
FindChirpSP
TaylorF2RedSpin
Eccentricity
EOBNR
PadeF1
NR_hdf5
IMRPhenomA
AmpCorPPN
TaylorF1
TaylorT3
TaylorT4
TaylorF2
NumRelNinja2
IMRPhenomB
PhenSpinTaylorRD
SpinTaylorFrameless
BCV
TaylorT1
SpinTaylor
SpinTaylorT3
TaylorT2
LAL_PNORDER_NUM_ORDER
int XLALSimInspiralREAL8WaveTaper(REAL8Vector *signalvec, LALSimInspiralApplyTaper bookends)
char char * XLALStringDuplicate(const char *s)
int XLALStringCaseCompare(const char *s1, const char *s2)
char * XLALStringToken(char **s, const char *delim, int empty)
static const INT4 r
void XLALDestroyREAL8TimeSeries(REAL8TimeSeries *series)
#define XLAL_ERROR(...)
int XLALPrintError(const char *fmt,...) _LAL_GCC_PRINTF_FORMAT_(1
XLAL_SUCCESS
XLAL_ENAME
XLAL_EFUNC
XLAL_FAILURE
SkyPosition position
REAL4TimeVectorSeries * h
The inspiral waveform parameter structure containing information about the waveform to be generated.
Definition: LALInspiral.h:205
Approximant approximant
Post-Newtonain approximant to be used in generating the wave (input)
Definition: LALInspiral.h:208
INT4 gpsNanoSeconds
This structure stores the parameters for constructing a restricted post-Newtonian waveform.
REAL4 psi
polarization angle (radians)
SkyPosition position
location of source on sky
REAL4VectorSequence * data
REAL4 * data
REAL8Sequence * data
LALUnit sampleUnits
LIGOTimeGPS epoch
REAL8 * data
CHAR taper[LIGOMETA_INSPIRALTAPER_MAX]
CHAR waveform[LIGOMETA_WAVEFORM_MAX]
CHAR numrel_data[LIGOMETA_STRING_MAX]
double inclination
double distance
double f_min
double deltaT