25#include <lal/TimeSeries.h>
26#include <lal/LALConstants.h>
27#include <lal/LALStdlib.h>
28#include <lal/LALSimInspiral.h>
29#include <lal/LALSimInspiralPrecess.h>
30#include <lal/LALSimSphHarmMode.h>
41#define ROTATEZ(angle, vx, vy, vz)\
42 tmp1 = vx*cos(angle) - vy*sin(angle);\
43 tmp2 = vx*sin(angle) + vy*cos(angle);\
47#define ROTATEY(angle, vx, vy, vz)\
48 tmp1 = vx*cos(angle) + vz*sin(angle);\
49 tmp2 = - vx*sin(angle) + vz*cos(angle);\
56 if ( (fabs(creal(val1) - creal(val2)) >
EPSILON) || (fabs(cimag(val1) - cimag(val2)) >
EPSILON) )
71int main (
int argc,
char **argv)
88 const REAL8 eta = m1/(m1+m2)*m2/(m1+m2);
91 const REAL8 chi1 = 0.8;
92 const REAL8 chi2 = 0.3;
155 for (
UINT4 idx=0;idx<Ntest;idx++) {
159 iota = idxr*
LAL_PI/NtestR;
160 iota2 = idxr*
LAL_PI/NtestR*0.4;
161 psi = (NtestR-idxr-1.)*2.*
LAL_PI/NtestR;
162 psi2 = (NtestR-idxr-1.)*2.*
LAL_PI/NtestR*0.4;
163 v = (idxr+1.)/NtestR;
232 errCode +=
XLALSimInspiralSpinPNMode2m(&hStart,vts,psi0ts,LNhx0,LNhy0,LNhz0,e1x0,e1y0,e1z0,S1x0,S1y0,S1z0,S2x0,S2y0,S2z0,m1,m2,dist,ampO);
233 errCode +=
XLALSimInspiralSpinPNMode2m(&hFin,vts,psits,LNhx,LNhy,LNhz,e1x,e1y,e1z,S1x,S1y,S1z,S2x,S2y,S2z,m1,m2,dist,ampO);
236 errCode +=
XLALSimInspiralSpinPNMode3m(&hStart,vts,psi0ts,LNhx0,LNhy0,LNhz0,e1x0,e1y0,e1z0,S1x0,S1y0,S1z0,S2x0,S2y0,S2z0,m1,m2,dist,ampO);
237 errCode +=
XLALSimInspiralSpinPNMode3m(&hFin,vts,psits,LNhx,LNhy,LNhz,e1x,e1y,e1z,S1x,S1y,S1z,S2x,S2y,S2z,m1,m2,dist,ampO);
240 errCode +=
XLALSimInspiralSpinPNMode4m(&hStart,vts,psi0ts,LNhx0,LNhy0,LNhz0,e1x0,e1y0,e1z0,S1x0,S1y0,S1z0,S2x0,S2y0,S2z0,m1,m2,dist,ampO);
241 errCode +=
XLALSimInspiralSpinPNMode4m(&hFin,vts,psits,LNhx,LNhy,LNhz,e1x,e1y,e1z,S1x,S1y,S1z,S2x,S2y,S2z,m1,m2,dist,ampO);
265 for(
UINT4 idx=0;idx<Ntest;idx++) {
289 for(
UINT4 idx=0;idx<Ntest;idx++) {
300 if ( (ret == 0) && (errCode == 0) ) {
301 fprintf(stdout,
"\n Precessing modes test passed.\n");
304 fprintf(stderr,
"\nFAILURE: %u Precessing modes test failed.\n", ret+errCode);
307 return ret + errCode ;
INT4 XLALSimInspiralSpinPNMode4m(SphHarmTimeSeries **hlm, REAL8TimeSeries *V, REAL8TimeSeries *Phi, REAL8TimeSeries *LNhx, REAL8TimeSeries *LNhy, REAL8TimeSeries *LNhz, REAL8TimeSeries *e1x, REAL8TimeSeries *e1y, REAL8TimeSeries *e1z, REAL8TimeSeries *S1x, REAL8TimeSeries *S1y, REAL8TimeSeries *S1z, REAL8TimeSeries *S2x, REAL8TimeSeries *S2y, REAL8TimeSeries *S2z, REAL8 m1, REAL8 m2, REAL8 distance, int ampO)
#define ROTATEZ(angle, vx, vy, vz)
int main(int argc, char **argv)
static int c_compare(COMPLEX16 val1, COMPLEX16 val2)
#define ROTATEY(angle, vx, vy, vz)
static int compare(REAL8 val1, REAL8 val2)
INT4 XLALSimInspiralSpinPNMode2m(SphHarmTimeSeries **hlm, REAL8TimeSeries *V, REAL8TimeSeries *Phi, REAL8TimeSeries *LNhx, REAL8TimeSeries *LNhy, REAL8TimeSeries *LNhz, REAL8TimeSeries *e1x, REAL8TimeSeries *e1y, REAL8TimeSeries *e1z, REAL8TimeSeries *S1x, REAL8TimeSeries *S1y, REAL8TimeSeries *S1z, REAL8TimeSeries *S2x, REAL8TimeSeries *S2y, REAL8TimeSeries *S2z, REAL8 m1, REAL8 m2, REAL8 distance, int ampO)
Computes the 5 l=2 modes of spherical harmonic decomposition of the post-Newtonian inspiral waveform ...
INT4 XLALSimInspiralSpinPNMode3m(SphHarmTimeSeries **hlm, REAL8TimeSeries *V, REAL8TimeSeries *Phi, REAL8TimeSeries *LNhx, REAL8TimeSeries *LNhy, REAL8TimeSeries *LNhz, REAL8TimeSeries *e1x, REAL8TimeSeries *e1y, REAL8TimeSeries *e1z, REAL8TimeSeries *S1x, REAL8TimeSeries *S1y, REAL8TimeSeries *S1z, REAL8TimeSeries *S2x, REAL8TimeSeries *S2y, REAL8TimeSeries *S2z, REAL8 m1, REAL8 m2, REAL8 distance, int ampO)
Computes all l=3 modes of spherical harmonic decomposition of the post-Newtonian inspiral waveform fo...
int XLALSimInspiralPrecessionRotateModesOut(SphHarmTimeSeries **hlm_out, SphHarmTimeSeries *hlm_in, const REAL8TimeSeries *alpha, const REAL8TimeSeries *beta, const REAL8TimeSeries *gam)
Takes in the h_lm spherical harmonic decomposed modes and rotates the modes by Euler angles alpha,...
int XLALSimAddModeAngleTimeSeries(REAL8TimeSeries *hplus, REAL8TimeSeries *hcross, COMPLEX16TimeSeries *hmode, REAL8TimeSeries *theta, REAL8TimeSeries *phi, int l, int m, int sym)
Multiplies a mode h(l,m) by a spin-2 weighted spherical harmonic to obtain hplus - i hcross,...
int XLALSimAddMode(REAL8TimeSeries *hplus, REAL8TimeSeries *hcross, COMPLEX16TimeSeries *hmode, REAL8 theta, REAL8 phi, int l, int m, int sym)
Multiplies a mode h(l,m) by a spin-2 weighted spherical harmonic to obtain hplus - i hcross,...
COMPLEX16TimeSeries * XLALSphHarmTimeSeriesGetMode(SphHarmTimeSeries *ts, UINT4 l, INT4 m)
Get the time series of a waveform's (l,m) spherical harmonic mode from a SphHarmTimeSeries linked lis...
void XLALDestroyCOMPLEX16TimeSeries(COMPLEX16TimeSeries *series)
REAL8TimeSeries * XLALCreateREAL8TimeSeries(const CHAR *name, const LIGOTimeGPS *epoch, REAL8 f0, REAL8 deltaT, const LALUnit *sampleUnits, size_t length)
const LALUnit lalStrainUnit
const LALUnit lalDimensionlessUnit
Structure to carry a collection of spherical harmonic modes in COMPLEX16 time series.