Loading [MathJax]/extensions/TeX/AMSsymbols.js
LALSimulation 6.2.0.1-8a6b96f
All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Modules Pages
PrecessWaveformIMRPhenomBTest.c
Go to the documentation of this file.
1/*
2 * Copyright (C) 2012 Chris Pankow, Evan Ochsner
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 Chris Pankow
22 *
23 * \file
24 *
25 * \brief Testing constant precession code on IMRPhenomB waveform.
26 */
27
28
29#include <lal/LALSimInspiralPrecess.h>
30
31int main(void){
32
33 FILE* h_ref = fopen("h_ref_PhenomB.txt", "w");
34 FILE* h_rot = fopen("h_rot_PhenomB.txt", "w");
35
36 int ret;
37 unsigned int i;
38
39 // Waveform parameters
40 double m1 = 2.0*LAL_MSUN_SI, m2 = 5.0*LAL_MSUN_SI;
41 double s1x = 0.0, s1y = 0.0, s1z = 0.0;
42 double s2x = 0.0, s2y = 0.0, s2z = 0.0;
43 double f_min = 40.0, f_ref = 0.0, dist = 1e6*LAL_PC_SI, inc = LAL_PI_4;
44 double lambda1 = 0.0, lambda2 = 0.0;
45 double phi = 0.0, dt = 1/16384.0;
46 LALSimInspiralWaveformFlags *waveFlags = NULL;
47 LALSimInspiralTestGRParam *nonGRparams = NULL;
48 int amplitudeOrder = 0, phaseOrder = 7;
50
51 double view_th = 0.0, view_ph = 0.0;
52 double Y_2_m2 = XLALSpinWeightedSphericalHarmonic( view_th + LAL_PI, view_ph, -2, 2, -2), Y_22 = XLALSpinWeightedSphericalHarmonic( view_th, view_ph, -2, 2, 2);
53
54 REAL8TimeSeries *hp, *hx;
55 hp = NULL; hx = NULL;
56
58 &hp, &hx,
59 m1, m2,
60 s1x, s1y, s1z,
61 s2x, s2y, s2z,
62 dist, inc,
63 phi, 0., 0., 0.,
64 dt, f_min, f_ref,
65 lambda1, lambda2, 0., 0.,
66 waveFlags,
67 nonGRparams,
68 amplitudeOrder, phaseOrder,
70 );
71 if( ret != XLAL_SUCCESS )
73
74 COMPLEX16TimeSeries *h_22, *h_2_2;
75 h_22 = NULL; h_2_2 = NULL;
76 // Initialize the h_lm REAL8 vectors
78 "h_{2-2}",
79 &(hp->epoch),
80 hp->f0,
81 hp->deltaT,
82 &(hp->sampleUnits),
83 hp->data->length
84 );
85
86 // Define h_{2-2}
87 for( i=0; i< hp->data->length; i++ ){
88 h_2_2->data->data[i] = (hp->data->data[i] + hx->data->data[i] * I) / Y_2_m2;
89 }
90
93 hp = NULL; hx = NULL;
94
95 inc=0; // +z axis
97 &hp, &hx,
98 m1, m2,
99 s1x, s1y, s1z,
100 s2x, s2y, s2z,
101 dist, inc,
102 phi, 0., 0., 0.,
103 dt, f_min, f_ref,
104 lambda1, lambda2, 0., 0.,
105 waveFlags,
106 NULL, // non-GR params
107 amplitudeOrder, phaseOrder,
109 );
110 if( ret != XLAL_SUCCESS )
112
113
114 // Initialize the h_lm REAL8 vectors
116 "h_{22}",
117 &(hp->epoch),
118 hp->f0,
119 hp->deltaT,
120 &(hp->sampleUnits),
121 hp->data->length
122 );
123
124 // Define h_{22}
125 for( i=0; i< hp->data->length; i++ ){
126 h_22->data->data[i] = (hp->data->data[i] + hx->data->data[i] * I) / Y_22;
127 }
128
129 // Write out reference waveform
130 REAL8 t0 = XLALGPSGetREAL8(&(hp->epoch));
131 for(i=0; i<hp->data->length; i++)
132 fprintf( h_ref, "%g %g %g\n", t0 + i * hp->deltaT,
133 hp->data->data[i], hx->data->data[i] );
134
136
137 SphHarmTimeSeries *ts = NULL;
138 ts = XLALSphHarmTimeSeriesAddMode( ts, h_22, 2, 2 );
139 ts = XLALSphHarmTimeSeriesAddMode( ts, h_2_2, 2, -2 );
140
142 &hp, &hx,
143 ts,
144 10, LAL_PI/4, 0,
145 0, LAL_PI/4 );
146 if( ret != XLAL_SUCCESS )
148
152
153 // Write out rotated waveform
154 for(i=0; i<hp->data->length; i++)
155 fprintf( h_rot, "%g %g %g\n", t0 + i * hp->deltaT,
156 hp->data->data[i], hx->data->data[i] );
157
158 // We're done.
161
162 fclose(h_ref);
163 fclose(h_rot);
164
165 return 0;
166}
int XLALSimInspiralChooseTDWaveformOLD(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, const REAL8 lambda1, const REAL8 lambda2, const REAL8 dQuadParam1, const REAL8 dQuadParam2, LALSimInspiralWaveformFlags *waveFlags, LALSimInspiralTestGRParam *nonGRparams, int amplitudeO, const int phaseO, const Approximant approximant)
int main(void)
#define fprintf
double dt
Definition: bh_ringdown.c:113
double i
Definition: bh_ringdown.c:118
#define LAL_MSUN_SI
#define LAL_PI
#define LAL_PC_SI
#define LAL_PI_4
double REAL8
Approximant
Enum that specifies the PN approximant to be used in computing the waveform.
@ IMRPhenomB
Time domain (non-precessing spins) inspiral-merger-ringdown waveforms generated from the inverse FFT ...
int XLALSimInspiralConstantPrecessionConeWaveform(REAL8TimeSeries **hp, REAL8TimeSeries **hx, SphHarmTimeSeries *h_lm, double precess_freq, double a, double phi_precess, double alpha_0, double beta_0)
Takes in the spherical harmonic decomposed modes as a strain time series and imposes the effect of a ...
void XLALSimInspiralDestroyWaveformFlags(LALSimInspiralWaveformFlags *waveFlags)
Destroy a LALSimInspiralWaveformFlags struct.
SphHarmTimeSeries * XLALSphHarmTimeSeriesAddMode(SphHarmTimeSeries *appended, const COMPLEX16TimeSeries *inmode, UINT4 l, INT4 m)
Prepend a node to a linked list of SphHarmTimeSeries, or create a new head.
void XLALDestroySphHarmTimeSeries(SphHarmTimeSeries *ts)
Delete list from current pointer to the end of the list.
COMPLEX16 XLALSpinWeightedSphericalHarmonic(REAL8 theta, REAL8 phi, int s, int l, int m)
void XLALDestroyCOMPLEX16TimeSeries(COMPLEX16TimeSeries *series)
COMPLEX16TimeSeries * XLALCreateCOMPLEX16TimeSeries(const CHAR *name, const LIGOTimeGPS *epoch, REAL8 f0, REAL8 deltaT, const LALUnit *sampleUnits, size_t length)
void XLALDestroyREAL8TimeSeries(REAL8TimeSeries *series)
#define XLAL_ERROR(...)
XLAL_SUCCESS
XLAL_EFUNC
REAL8 XLALGPSGetREAL8(const LIGOTimeGPS *epoch)
COMPLEX16Sequence * data
COMPLEX16 * data
Linked list of any number of parameters for testing GR.
REAL8Sequence * data
LALUnit sampleUnits
LIGOTimeGPS epoch
REAL8 * data
Structure to carry a collection of spherical harmonic modes in COMPLEX16 time series.
double f_min
Definition: unicorn.c:22