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
LALInspiralWave3.c
Go to the documentation of this file.
1/*
2* Copyright (C) 2007 Stas Babak, David Churches, B.S. Sathyaprakash, Craig Robinson , Thomas Cokelaer
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 Sathyaprakash, B. S.
22 * \file
23 * \ingroup LALInspiral_h
24 *
25 * \brief These modules generate a time-domain chirp waveform of type #TaylorT3.
26 *
27 * ### Prototypes ###
28 *
29 * <tt>LALInspiralWave3()</tt>
30 * <ul>
31 * <li> \c output: Output containing the inspiral waveform.</li>
32 * <li> \c params: Input containing binary chirp parameters.</li>
33 * </ul>
34 *
35 * <tt>XLALInspiralWave3Templates()</tt>
36 * <ul>
37 * <li> \c output1: Output containing the 0-phase inspiral waveform.</li>
38 * <li> \c output2: Output containing the \f$\pi/2\f$-phase inspiral waveform.</li>
39 * <li> \c params: Input containing binary chirp parameters.</li>
40 * </ul>
41 *
42 * ### Description ###
43 *
44 * LALInspiralWave3() generates #TaylorT3 approximant which
45 * corresponds to the case wherein
46 * the phase of the waveform is given as an explicit function of time
47 * as in \eqref{eq_InspiralWavePhase3}.
48 *
49 * XLALInspiralWave3Templates() simultaneously generates
50 * two inspiral waveforms and the two differ in
51 * phase by \f$\pi/2\f$.
52 *
53 * ### Algorithm ###
54 *
55 *
56 * ### Uses ###
57 *
58 * \code
59 * LALInspiralParameterCalc()
60 * LALInspiralChooseModel()
61 * LALInspiralSetup()
62 * LALInspiralPhasing3 (via expnFunc)()
63 * LALInspiralFrequency3 (via expnFunc)()
64 * \endcode
65 *
66 */
67
68#include <math.h>
69#include <lal/LALStdlib.h>
70#include <lal/LALInspiral.h>
71#include <lal/FindRoot.h>
72#include <lal/Units.h>
73#include <lal/SeqFactories.h>
74
75
76typedef struct
77{
78 REAL8 (*func)(REAL8 tC, expnCoeffs *ak);
80}
82
83static REAL8 XLALInspiralFrequency3Wrapper(REAL8 tC, void *pars);
84
85static int
87 REAL4Vector *output1,
88 REAL4Vector *output2,
89 REAL4Vector *h,
91 REAL4Vector *ff,
92 REAL8Vector *phi,
94 InspiralInit *paramsInit
95 );
96
97void
100 REAL4Vector *output,
102 )
103{
104 XLAL_PRINT_DEPRECATION_WARNING("XLALInspiralWave3");
107
110
112 RETURN(status);
113}
114
115int
117 REAL4Vector *output,
119 )
120{
121 INT4 count;
122 InspiralInit paramsInit;
123
124 if (output == NULL)
126 if (output->data == NULL)
128 if (params == NULL)
130 if (params->nStartPad < 0)
132 if (params->fLower <= 0)
134 if (params->tSampling <= 0)
136
139 if (XLALInspiralSetup(&(paramsInit.ak), params))
141 if (XLALInspiralChooseModel(&(paramsInit.func), &(paramsInit.ak), params))
143
144 if (params->totalMass <= 0.)
146 if (params->eta < 0.)
148
149 memset( output->data, 0, output->length * sizeof(REAL4) );
150
151 /* Call the engine function */
152 count = XLALInspiralWave3Engine(output, NULL, NULL,
153 NULL, NULL, NULL, params, &paramsInit);
154 if (count < 0)
156
157 return XLAL_SUCCESS;
158}
159
161{
163 REAL8 freq, f;
164
165 in = (ChirptimeFromFreqIn *) pars;
166 freq = in->func(tC, &(in->ak));
167 if (XLAL_IS_REAL8_FAIL_NAN(freq))
169 f = freq - in->ak.f0;
170
171 /*
172 fprintf(stderr, "Here freq=%e f=%e tc=%e f0=%e\n", freq, *f, tC, in->ak.f0);
173 */
174
175 return f;
176}
177
178int
180 REAL4Vector *output1,
181 REAL4Vector *output2,
183 )
184{
185 INT4 count;
186
187 InspiralInit paramsInit;
188
189 if (output1 == NULL)
191 if (output2 == NULL)
193 if (output1->data == NULL)
195 if (output2->data == NULL)
197 if (params == NULL)
199 if (params->nStartPad < 0)
201 if (params->fLower <= 0)
203 if (params->tSampling <= 0)
205
208 if (XLALInspiralSetup(&(paramsInit.ak), params))
210 if (XLALInspiralChooseModel(&(paramsInit.func), &(paramsInit.ak), params))
212
213 if (params->totalMass <= 0.)
215 if (params->eta < 0.)
217
218 /* Initialise the waveforms to zero */
219 memset(output1->data, 0, output1->length * sizeof(REAL4));
220 memset(output2->data, 0, output2->length * sizeof(REAL4));
221
222 /* Call the engine function */
223 count = XLALInspiralWave3Engine(output1, output2, NULL,
224 NULL, NULL, NULL, params, &paramsInit);
225 if (count < 0)
227
228 return XLAL_SUCCESS;
229}
230
231int
233 CoherentGW *waveform,
235 PPNParamStruc *ppnParams
236 )
237{
238 INT4 count;
239 UINT4 i;
240 REAL4Vector *h=NULL;
241 REAL4Vector *a=NULL;
242 REAL4Vector *ff=NULL ;
243 REAL8Vector *phiv=NULL;
244
245 REAL8 phiC;/* phase at coalescence */
246 CHAR message[256];
247
248
249 InspiralInit paramsInit;
250
251 /** -- -- */
252
253 /* Make sure parameter and waveform structures exist. */
254 if (params == NULL)
256 if (waveform == NULL)
258 if (waveform->h != NULL)
260 if (waveform->a != NULL)
262 if (waveform->f != NULL)
264 if (waveform->phi != NULL)
266 if (waveform->shift != NULL)
268
269 params->ampOrder = (LALPNOrder) 0;
270 sprintf(message, "WARNING: Amp Order has been reset to %d", params->ampOrder);
271 XLALPrintInfo("%s", message);
272 /* Compute some parameters*/
273 if (XLALInspiralInit(params, &paramsInit))
275
276 if (paramsInit.nbins==0)
277 {
278 /* FIXME: is this the correct thing to return? */
279 return XLAL_SUCCESS;
280 }
281
282 /* Now we can allocate memory and vector for coherentGW structure*/
283 ff = XLALCreateREAL4Vector(paramsInit.nbins);
284 if (ff == NULL)
286 a = XLALCreateREAL4Vector(2*paramsInit.nbins);
287 if (a == NULL)
289 phiv = XLALCreateREAL8Vector(paramsInit.nbins);
290 if (phiv == NULL)
292
293 /* By default the waveform is empty */
294
295 memset(ff->data, 0, ff->length * sizeof(REAL4));
296 memset(a->data, 0, a->length * sizeof(REAL4));
297 memset(phiv->data, 0, phiv->length * sizeof(REAL8));
298
299 if( params->ampOrder )
300 {
301 h = XLALCreateREAL4Vector(2*paramsInit.nbins);
302 if (h == NULL)
304 memset(h->data, 0, h->length * sizeof(REAL4));
305 }
306
307 /* Call the engine function */
308 count = XLALInspiralWave3Engine(NULL, NULL, h, a, ff, phiv, params, &paramsInit);
309 if (count < 0)
310 {
314 if( h )
315 {
317 }
319 }
320
321 /* Check an empty waveform hasn't been returned */
322 for (i = 0; i < phiv->length; i++)
323 {
324 if (phiv->data[i] != 0.0) break;
325 /* If the waveform returned is empty, return now */
326 if (i == phiv->length - 1)
327 {
331 if( h )
332 {
334 }
335
336 /* FIXME: is this the correct thing to return? */
337 return XLAL_SUCCESS;
338 }
339 }
340
341 /* if ( (phase/2./LAL_PI) < 2. ){
342 sprintf(message, "The waveform has only %lf cycles; we don't keep waveform with less than 2 cycles.",
343 (double)phase/2./(double)LAL_PI );
344 LALWarning(status, message);
345
346
347 }
348 else*/
349 {
350
351 /*wrap the phase vector*/
352 phiC = phiv->data[count-1] ;
353 for (i=0; i<(UINT4)count;i++)
354 {
355 phiv->data[i] = phiv->data[i] -phiC + ppnParams->phi;
356 }
357 /* Allocate the waveform structures. */
358 waveform->a = (REAL4TimeVectorSeries *) XLALMalloc( sizeof(REAL4TimeVectorSeries) );
359 if ( waveform->a == NULL )
361 memset( waveform->a, 0, sizeof(REAL4TimeVectorSeries) );
362
363 waveform->f = (REAL4TimeSeries *) LALMalloc( sizeof(REAL4TimeSeries) );
364 if ( waveform->f == NULL )
365 {
366 XLALFree( waveform->a );
367 waveform->a = NULL;
369 }
370 memset( waveform->f, 0, sizeof(REAL4TimeSeries) );
371
372 waveform->phi = (REAL8TimeSeries *) LALMalloc( sizeof(REAL8TimeSeries) );
373 if ( waveform->phi == NULL )
374 {
375 XLALFree( waveform->a );
376 waveform->a = NULL;
377 XLALFree( waveform->f );
378 waveform->f = NULL;
380 }
381 memset( waveform->phi, 0, sizeof(REAL8TimeSeries) );
382
383 waveform->a->data = XLALCreateREAL4VectorSequence(count, 2);
384 if (waveform->a->data == NULL)
386 waveform->f->data = XLALCreateREAL4Vector(count);
387 if (waveform->f->data == NULL)
389 waveform->phi->data = XLALCreateREAL8Vector(count);
390 if (waveform->phi->data == NULL)
392
393 memcpy(waveform->f->data->data , ff->data, count*(sizeof(REAL4)));
394 memcpy(waveform->a->data->data , a->data, 2*count*(sizeof(REAL4)));
395 memcpy(waveform->phi->data->data ,phiv->data, count*(sizeof(REAL8)));
396
397 waveform->a->deltaT = waveform->f->deltaT = waveform->phi->deltaT = 1./params->tSampling;
398
399 waveform->a->sampleUnits = lalStrainUnit;
400 waveform->f->sampleUnits = lalHertzUnit;
402 waveform->position = ppnParams->position;
403 waveform->psi = ppnParams->psi;
404
405 snprintf( waveform->a->name, LALNameLength, "T3 inspiral amplitudes" );
406 snprintf( waveform->f->name, LALNameLength, "T3 inspiral frequency" );
407 snprintf( waveform->phi->name, LALNameLength, "T3 inspiral phase" );
408
409 /* --- fill some output ---*/
410
411 ppnParams->tc = (double)(count-1) / params->tSampling ;
412 ppnParams->length = count;
413 ppnParams->dfdt = ((REAL4)(waveform->f->data->data[count-1] - waveform->f->data->data[count-2])) * ppnParams->deltaT;
414 ppnParams->fStop = params->fFinal;
416 ppnParams->termDescription = GENERATEPPNINSPIRALH_MSGEFSTOP;
417
418 ppnParams->fStart = ppnParams->fStartIn;
419
420 if( params->ampOrder )
421 {
422 waveform->h = (REAL4TimeVectorSeries *) XLALMalloc( sizeof(REAL4TimeVectorSeries) );
423 if ( waveform->h == NULL )
425 memset( waveform->h, 0, sizeof(REAL4TimeVectorSeries) );
426
427 waveform->h->data = XLALCreateREAL4VectorSequence(count, 2);
428 if ( waveform->h->data == NULL )
430 memcpy(waveform->h->data->data , h->data, 2*count*(sizeof(REAL4)));
431 waveform->h->deltaT = 1./params->tSampling;
432 waveform->h->sampleUnits = lalStrainUnit;
433 snprintf( waveform->h->name, LALNameLength, "T3 inspiral polarizations" );
435 h = NULL;
436 }
437 } /*end of coherentGW storage */
438
439
440 /* --- free memory --- */
444
445 return XLAL_SUCCESS;
446}
447
448
449/* Engine function used to generate the waveforms */
450static int
452 REAL4Vector *output1,
453 REAL4Vector *output2,
454 REAL4Vector *h,
455 REAL4Vector *a,
456 REAL4Vector *ff,
457 REAL8Vector *phiv,
459 InspiralInit *paramsInit
460 )
461
462{
463 INT4 i, startShift, count;
464 REAL8 dt, fu, eta, tc, totalMass, t, td, c1, phi0, phi1, phi;
465 REAL8 v, v2, f, fHigh, tmax, fOld, phase, omega;
466 REAL8 xmin,xmax,xacc;
467 REAL8 (*frequencyFunction)(REAL8, void *);
468 expnFunc func;
469 expnCoeffs ak;
470 ChirptimeFromFreqIn timeIn;
471 void *pars;
472
473 REAL8 temp, tempMax=0, tempMin = 0;
474
475 /* Only used in injection case */
476 REAL8 unitHz = 0.;
477 REAL8 cosI = 0.;/* cosine of system inclination */
478 REAL8 apFac = 0., acFac = 0.;/* extra factor in plus and cross amplitudes */
479
480
481 ak = paramsInit->ak;
482 func = paramsInit->func;
483
484 if (output2 || a)
485 params->nStartPad = 0; /* must be zero for templates and injections */
486
487 eta = params->eta; /* Symmetric mass ratio */
488 /* Only in injection case.. */
489 unitHz = params->totalMass*LAL_MTSUN_SI*(REAL8)LAL_PI;
490 cosI = cos( params->inclination );
491 apFac = acFac = -2.0 * eta *params->totalMass * LAL_MRSUN_SI/params->distance;
492 apFac *= 1.0 + cosI*cosI;
493 acFac *= 2.0 * cosI;
494
495 dt = 1.0/(params->tSampling); /* sampling rate */
496 fu = params->fCutoff; /* upper frequency cutoff */
497 phi = params->startPhase; /* initial phase */
498 startShift = params->nStartPad; /* number of zeros at the start of the wave */
499
500
501 tc=params->tC; /* Instant of coalescence of the compact objects */
502 totalMass = (params->totalMass)*LAL_MTSUN_SI; /* mass of the system in seconds */
503
504
505/*
506 If flso is less than the user inputted upper frequency cutoff fu,
507 then the waveforn is truncated at f=flso.
508*/
509
510 if (fu)
511 fHigh = (fu < ak.flso) ? fu : ak.flso;
512 else
513 fHigh = ak.flso;
514
515/*
516 Check that the highest frequency is less than half the sampling frequency -
517 the Nyquist theorum
518*/
519
520 if (fHigh >= 0.5/dt)
521 {
522 XLALPrintError("fHigh must be less than Nyquist frequency\n");
524 }
525
526 if (fHigh <= params->fLower)
527 {
528 XLALPrintError("fHigh must be larger than fLower\n");
530 }
531
532/* Here's the part which calculates the waveform */
533
534 c1 = eta/(5.*totalMass);
535
536 i = startShift;
537
538 /*
539 * In Jan 2003 we realized that the tC determined as a sum of chirp times is
540 * not quite the tC that should enter the definition of Theta in the expression
541 * for the frequency as a function of time (see DIS3, 2000). This is because
542 * chirp times are obtained by inverting t(f). Rather tC should be obtained by
543 * solving the equation f0 - f(tC) = 0. This is what is implemented below.
544 */
545
546 timeIn.func = func.frequency3;
547 timeIn.ak = ak;
548 frequencyFunction = &XLALInspiralFrequency3Wrapper;
549 xmin = c1*params->tC/2.;
550 xmax = c1*params->tC*2.;
551 xacc = 1.e-6;
552 pars = (void*) &timeIn;
553 /* tc is the instant of coalescence */
554
555 xmax = c1*params->tC*3 + 5.; /* we add 5 so that if tC is small then xmax
556 is always greater than a given value (here 5)*/
557
558 /* for x in [xmin, xmax], we search the value which gives the max frequency.
559 and keep the corresponding rootIn.xmin. */
560
561
562
563 for (tc = c1*params->tC/1000.; tc < xmax; tc+=c1*params->tC/1000.){
564 temp = XLALInspiralFrequency3Wrapper(tc , pars);
565 if (XLAL_IS_REAL8_FAIL_NAN(temp))
567 if (temp > tempMax) {
568 xmin = tc;
569 tempMax = temp;
570 }
571 if (temp < tempMin) {
572 tempMin = temp;
573 }
574 }
575
576 /* if we have found a value positive then everything should be fine in the
577 BissectionFindRoot function */
578 if (tempMax > 0 && tempMin < 0){
579 tc = XLALDBisectionFindRoot (frequencyFunction, xmin, xmax, xacc, pars);
582 }
583 else if (a)
584 {
585 /* Otherwise we return an empty waveform for injection */
586 return 0;
587 }
588 else
589 {
590 /* Or abort if not injection */
591 XLALPrintError("Can't find good bracket for BisectionFindRoot\n");
593 }
594
595 tc /= c1;
596
597 tc += params->startTime; /* Add user given startTime to instant of
598 coalescence of the compact objects */
599
600 t=0.0;
601 td = c1*(tc-t);
602 phase = func.phasing3(td, &ak);
603 if (XLAL_IS_REAL8_FAIL_NAN(phase))
605 f = func.frequency3(td, &ak);
608 phi0=-phase+phi;
609 phi1=phi0+LAL_PI_2;
610
611 count = 0;
612 tmax = tc - dt;
613 fOld = 0.0;
614
615/* We stop if any of the following conditions fail */
616
617 while (f < fHigh && t < tmax && f > fOld)
618 {
619 /* Check we don't write past the end of the vector */
620 if ((output1 && ((UINT4)i >= output1->length)) || (ff && ((UINT4)count >= ff->length)))
621 {
622 XLALPrintError("Attempting to write beyond the end of vector\n");
624 }
625
626 fOld = f;
627 v = pow(f*LAL_PI*totalMass, (1./3.));
628 v2 = v*v;
629 if (output1)
630 {
631
632
633 /*
634 output1->data[i] = LALInspiralHPlusPolarization( phase+phi0, v, params );
635 */
636 output1->data[i] = (REAL4) (apFac * v2) * cos(phase+phi0);
637 if (output2)
638 /*
639 output2->data[i] = LALInspiralHCrossPolarization( phase+phi1, v, params );
640 */
641 output2->data[i] = (REAL4) (apFac * v2) * cos(phase+phi1);
642 }
643 else if (a)
644 {
645 int ice, ico;
646 ice = 2*count;
647 ico = ice + 1;
648 omega = v*v*v;
649
650 ff->data[count] = (REAL4)(omega/unitHz);
651 a->data[ice] = (REAL4)(apFac * v2 );
652 a->data[ico] = (REAL4)(acFac * v2 );
653 phiv->data[count] = (REAL8)(phase );
654
655 if (h)
656 {
657 h->data[ice] = LALInspiralHPlusPolarization( phase, v, params );
658 h->data[ico] = LALInspiralHCrossPolarization( phase, v, params );
659 }
660 }
661 ++i;
662 ++count;
663 t=count*dt;
664 td = c1*(tc-t);
665 phase = func.phasing3(td, &ak);
666 if (XLAL_IS_REAL8_FAIL_NAN(phase))
668 f = func.frequency3(td, &ak);
671 }
672 params->fFinal = fOld;
673 if (output1 && !output2) params->tC = t;
674/*
675 fprintf(stderr, "%e %e\n", f, fHigh);
676*/
677
678 return count;
679}
int XLALInspiralParameterCalc(InspiralTemplate *params)
REAL4 LALInspiralHCrossPolarization(REAL8 phase, REAL8 v, InspiralTemplate *params)
int XLALInspiralSetup(expnCoeffs *ak, InspiralTemplate *params)
int XLALInspiralInit(InspiralTemplate *params, InspiralInit *paramsInit)
int XLALInspiralChooseModel(expnFunc *func, expnCoeffs *ak, InspiralTemplate *params)
REAL4 LALInspiralHPlusPolarization(REAL8 phase, REAL8 v, InspiralTemplate *params)
static int XLALInspiralWave3Engine(REAL4Vector *output1, REAL4Vector *output2, REAL4Vector *h, REAL4Vector *a, REAL4Vector *ff, REAL8Vector *phi, InspiralTemplate *params, InspiralInit *paramsInit)
static REAL8 XLALInspiralFrequency3Wrapper(REAL8 tC, void *pars)
int XLALInspiralWave3Templates(REAL4Vector *output1, REAL4Vector *output2, InspiralTemplate *params)
void LALInspiralWave3(LALStatus *status, REAL4Vector *output, InspiralTemplate *params)
int XLALInspiralWave3(REAL4Vector *output, InspiralTemplate *params)
int XLALInspiralWave3ForInjection(CoherentGW *waveform, InspiralTemplate *params, PPNParamStruc *ppnParams)
#define LALMalloc(n)
const double c1
#define ATTATCHSTATUSPTR(statusptr)
#define DETATCHSTATUSPTR(statusptr)
#define INITSTATUS(statusptr)
#define RETURN(statusptr)
#define ABORTXLAL(sp)
double i
REAL4VectorSequence * XLALCreateREAL4VectorSequence(UINT4 length, UINT4 veclen)
REAL8 XLALDBisectionFindRoot(REAL8(*y)(REAL8, void *), REAL8 xmin, REAL8 xmax, REAL8 xacc, void *params)
#define GENERATEPPNINSPIRALH_EFSTOP
Reached requested termination frequency.
#define LAL_PI_2
#define LAL_PI
#define LAL_MTSUN_SI
#define LAL_MRSUN_SI
double REAL8
char CHAR
uint32_t UINT4
int32_t INT4
float REAL4
LALNameLength
void * XLALMalloc(size_t n)
void XLALFree(void *p)
LALPNOrder
static const INT4 a
const LALUnit lalStrainUnit
const LALUnit lalHertzUnit
const LALUnit lalDimensionlessUnit
REAL4Vector * XLALCreateREAL4Vector(UINT4 length)
void XLALDestroyREAL4Vector(REAL4Vector *vector)
REAL8Vector * XLALCreateREAL8Vector(UINT4 length)
void XLALDestroyREAL8Vector(REAL8Vector *vector)
#define XLAL_ERROR_REAL8(...)
int int int XLALPrintInfo(const char *fmt,...) _LAL_GCC_PRINTF_FORMAT_(1
#define XLAL_ERROR(...)
int XLALPrintError(const char *fmt,...) _LAL_GCC_PRINTF_FORMAT_(1
#define XLAL_PRINT_DEPRECATION_WARNING(replacement)
#define XLAL_IS_REAL8_FAIL_NAN(val)
XLAL_EBADLEN
XLAL_ENOMEM
XLAL_SUCCESS
XLAL_EFAULT
XLAL_EFUNC
XLAL_EDOM
XLAL_EFAILED
REAL8(* func)(REAL8 tC, expnCoeffs *ak)
REAL4TimeSeries * shift
SkyPosition position
REAL4TimeVectorSeries * a
REAL8TimeSeries * phi
REAL4TimeVectorSeries * h
REAL4TimeSeries * f
expnFunc func
Definition: LALInspiral.h:680
expnCoeffs ak
Definition: LALInspiral.h:679
The inspiral waveform parameter structure containing information about the waveform to be generated.
Definition: LALInspiral.h:205
This structure stores the parameters for constructing a restricted post-Newtonian waveform.
const CHAR * termDescription
The termination code description (above)
INT4 termCode
The termination condition (above) that stopped computation of the waveform.
REAL4 dfdt
The maximum value of encountered over any timestep used in generating the waveform.
REAL4 phi
The phase at coalescence (or arbitrary reference phase for a post -Newtonian approximation),...
REAL4 psi
polarization angle (radians)
REAL4 fStart
The actual starting frequency of the waveform, in Hz (normally close but not identical to fStartIn)
SkyPosition position
location of source on sky
UINT4 length
The length of the generated waveform.
REAL4 fStop
The frequency at the termination of the waveform, in Hz.
REAL8 tc
The time from the start of the waveform to coalescence (in the point-mass approximation),...
REAL4 fStartIn
The requested starting frequency of the waveform, in Hz.
REAL8 deltaT
The requested sampling interval of the waveform, in s.
CHAR name[LALNameLength]
REAL4Sequence * data
LALUnit sampleUnits
CHAR name[LALNameLength]
REAL4VectorSequence * data
REAL4 * data
REAL8Sequence * data
LALUnit sampleUnits
CHAR name[LALNameLength]
REAL8 * data
This structure contains various post-Newtonian and P-approximant expansion coefficients; the meanings...
Definition: LALInspiral.h:399
REAL8 flso
Definition: LALInspiral.h:487
Structure to hold the pointers to the generic functions defined above.
Definition: LALInspiral.h:546
InspiralPhasing3 * phasing3
Definition: LALInspiral.h:551
InspiralFrequency3 * frequency3
Definition: LALInspiral.h:552
double inclination
double distance
char output[FILENAME_MAX]