LAL  6.21.0.1-1696c80
DetResponse.c
Go to the documentation of this file.
1 /*
2 * Copyright (C) 2007 David Chin, Jolien Creighton, Kipp Cannon, Peter Shawhan
3 * Copyright (C) 2012 Matthew Pitkin
4 *
5 * This program is free software; you can redistribute it and/or modify
6 * it under the terms of the GNU General Public License as published by
7 * the Free Software Foundation; either version 2 of the License, or
8 * (at your option) any later version.
9 *
10 * This program is distributed in the hope that it will be useful,
11 * but WITHOUT ANY WARRANTY; without even the implied warranty of
12 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13 * GNU General Public License for more details.
14 *
15 * You should have received a copy of the GNU General Public License
16 * along with with program; see the file COPYING. If not, write to the
17 * Free Software Foundation, Inc., 59 Temple Place, Suite 330, Boston,
18 * MA 02111-1307 USA
19 */
20 
21 #include <math.h>
22 #include <gsl/gsl_sf_trig.h>
23 #include <lal/LALStdlib.h>
24 #include <lal/LALStdio.h>
25 #include <lal/LALError.h>
26 #include <lal/Units.h>
27 #include <lal/SeqFactories.h>
28 #include <lal/Date.h>
29 #include <lal/SkyCoordinates.h>
30 #include <lal/DetResponse.h>
31 #include <lal/TimeSeries.h>
32 #include <lal/XLALError.h>
33 
34 /**
35  * An implementation of the detector response formulae in Anderson et al
36  * PRD 63 042003 (2001) \cite ABCF2001 .
37  *
38  * Computes F+ and Fx for a source at a specified sky position,
39  * polarization angle, and sidereal time. Also requires the detector's
40  * response matrix which is defined by Eq. (B6) of [ABCF] using either
41  * Table 1 of \cite ABCF2001 or Eqs. (B11)--(B17) to compute the arm
42  * direction unit vectors.
43  */
45  double *fplus, /**< Returned value of F+ */
46  double *fcross, /**< Returned value of Fx */
47  const REAL4 D[3][3], /**< Detector response 3x3 matrix */
48  const double ra, /**< Right ascention of source (radians) */
49  const double dec, /**< Declination of source (radians) */
50  const double psi, /**< Polarization angle of source (radians) */
51  const double gmst /**< Greenwich mean sidereal time (radians) */
52 )
53 {
54  int i;
55  double X[3];
56  double Y[3];
57 
58  /* Greenwich hour angle of source (radians). */
59  const double gha = gmst - ra;
60 
61  /* pre-compute trig functions */
62  const double cosgha = cos(gha);
63  const double singha = sin(gha);
64  const double cosdec = cos(dec);
65  const double sindec = sin(dec);
66  const double cospsi = cos(psi);
67  const double sinpsi = sin(psi);
68 
69  /* Eq. (B4) of [ABCF]. Note that dec = pi/2 - theta, and gha =
70  * -phi where theta and phi are the standard spherical coordinates
71  * used in that paper. */
72  X[0] = -cospsi * singha - sinpsi * cosgha * sindec;
73  X[1] = -cospsi * cosgha + sinpsi * singha * sindec;
74  X[2] = sinpsi * cosdec;
75 
76  /* Eq. (B5) of [ABCF]. Note that dec = pi/2 - theta, and gha =
77  * -phi where theta and phi are the standard spherical coordinates
78  * used in that paper. */
79  Y[0] = sinpsi * singha - cospsi * cosgha * sindec;
80  Y[1] = sinpsi * cosgha + cospsi * singha * sindec;
81  Y[2] = cospsi * cosdec;
82 
83  /* Now compute Eq. (B7) of [ABCF] for each polarization state, i.e.,
84  * with s+=1 and sx=0 to get F+, with s+=0 and sx=1 to get Fx */
85  *fplus = *fcross = 0.0;
86  for(i = 0; i < 3; i++) {
87  const double DX = D[i][0] * X[0] + D[i][1] * X[1] + D[i][2] * X[2];
88  const double DY = D[i][0] * Y[0] + D[i][1] * Y[1] + D[i][2] * Y[2];
89  *fplus += X[i] * DX - Y[i] * DY;
90  *fcross += X[i] * DY + Y[i] * DX;
91  }
92 }
93 
94 
95 /**
96  *
97  * An implementation of the detector response for all six tensor, vector and
98  * scalar polarisation modes of general metric theories of gravity. We follow
99  * the convention of \cite Blaut2012 (this is also equivalent to the
100  * Equations in \cite ABCF2001 and \cite Nishizawa2009 in albeit with a
101  * rotated set of coordinates), but with \cite Blaut2012's \f$\theta = \pi/2 -
102  * dec\f$ and \f$\phi = ra-GMST\f$ rather than the gha = gmst - ra used here.
103  *
104  * The values computed are the tensor mode response, F+ and Fx ("cross"), the
105  * scalar breathing and longitudinal modes, Fb and Fl, and the vector "x" and
106  * "y" modes, Fx ("x") and Fy, for a source at a specified sky position,
107  * polarization angle, and sidereal time. Also requires the detector's
108  * response matrix which is defined by Eq. (B6) of \cite ABCF2001 using either
109  * Table 1 of \cite ABCF2001 or Eqs. (B11)--(B17) to compute the arm
110  * direction unit vectors.
111  */
113  double *fplus, /**< Returned value of F+ */
114  double *fcross, /**< Returned value of Fx (cross) */
115  double *fb, /**< Returned value of Fb (breathing mode) */
116  double *fl, /**< Returned value of Fl (scalar longitudinal) */
117  double *fx, /**< Returned value of Fx ("x" vector mode) */
118  double *fy, /**< Returned value of Fy (y vector mode) */
119  const REAL4 D[3][3], /**< Detector response 3x3 matrix */
120  const double ra, /**< Right ascention of source (radians) */
121  const double dec, /**< Declination of source (radians) */
122  const double psi, /**< Polarization angle of source (radians) */
123  const double gmst /**< Greenwich mean sidereal time (radians) */
124 )
125 {
126  int i;
127  double X[3];
128  double Y[3];
129  double Z[3];
130 
131  /* Greenwich hour angle of source (radians). */
132  const double gha = gmst - ra;
133 
134  /* pre-compute trig functions */
135  const double cosgha = cos(gha);
136  const double singha = sin(gha);
137  const double cosdec = cos(dec);
138  const double sindec = sin(dec);
139  const double cospsi = cos(psi);
140  const double sinpsi = sin(psi);
141 
142  /* Eq. (B4) of [ABCF]. Note that dec = pi/2 - theta, and gha =
143  * -phi where theta and phi are the standard spherical coordinates
144  * used in that paper. */
145  X[0] = -cospsi * singha - sinpsi * cosgha * sindec;
146  X[1] = -cospsi * cosgha + sinpsi * singha * sindec;
147  X[2] = sinpsi * cosdec;
148 
149  /* Eq. (B5) of [ABCF]. Note that dec = pi/2 - theta, and gha =
150  * -phi where theta and phi are the standard spherical coordinates
151  * used in that paper. */
152  Y[0] = sinpsi * singha - cospsi * cosgha * sindec;
153  Y[1] = sinpsi * cosgha + cospsi * singha * sindec;
154  Y[2] = cospsi * cosdec;
155 
156  /* Eqns from [Blaut2012] - but converted from Blaut's theta = pi/2 - dec
157  * given cos(dec) = sin(theta) and cos(theta) = sin(dec), and Blaut's phi = ra
158  * - gmst, so cos(phi) = cos(gha) and sin(phi) = -sin(gha). This is consistent
159  * with the convention in XLALComputeDetAMResponse.
160  */
161  Z[0] = -cosdec * cosgha;
162  Z[1] = cosdec * singha;
163  Z[2] = -sindec;
164 
165  /* Now compute Eq. (B7) of [ABCF] for each polarization state, i.e.,
166  * with s+=1 and sx=0 to get F+, with s+=0 and sx=1 to get Fx */
167  *fplus = *fcross = *fb = *fl = *fx = *fy = 0.0;
168  for(i = 0; i < 3; i++) {
169  const double DX = D[i][0] * X[0] + D[i][1] * X[1] + D[i][2] * X[2];
170  const double DY = D[i][0] * Y[0] + D[i][1] * Y[1] + D[i][2] * Y[2];
171  const double DZ = D[i][0] * Z[0] + D[i][1] * Z[1] + D[i][2] * Z[2];
172 
173  *fplus += X[i] * DX - Y[i] * DY;
174  *fcross += X[i] * DY + Y[i] * DX;
175 
176  /* scalar and vector modes from [Nishizawa2009]
177  * also Eq. (13.96) in the textbook by Poisson & Will */
178  *fb += X[i] * DX + Y[i] * DY;
179  *fl += Z[i] * DZ;
180  *fx += X[i] * DZ + Z[i] * DX;
181  *fy += Y[i] * DZ + Z[i] * DY;
182  }
183 }
184 
185 
186 /*
187  *
188  * beta = pi f L / c
189  * mu = k . u
190  *
191  * @sa
192  * John T. Whelan, "Higher-Frequency Corrections to Stochastic Formulae",
193  * LIGO-T070172.
194  * @sa
195  * Louis J. Rubbo, Neil J. Cornish, and Olivier Poujade, "Forward modeling of
196  * space-borne gravitational wave detectors", Phys. Rev. D 69, 082003 (2004);
197  * arXiv:gr-qc/0311069.
198  * @sa
199  * Malik Rakhmanov, "Response of LIGO to Gravitational Waves at High
200  * Frequencies and in the Vicinity of the FSR (37.5 kHz)", LIGO-T060237.
201  */
203 {
204  COMPLEX16 ans;
205  double Pibeta = LAL_PI * beta;
206  ans = cexp(I * Pibeta * (1.0 - mu)) * gsl_sf_sinc(beta * (1.0 + mu));
207  ans += cexp(-I * Pibeta * (1.0 + mu)) * gsl_sf_sinc(beta * (1.0 - mu));
208  ans *= 0.5;
209  return ans;
210 }
211 
212 
213 static void getarm(double u[3], double alt, double azi, double lat, double lon)
214 {
215  double cosalt = cos(alt);
216  double sinalt = sin(alt);
217  double cosazi = cos(azi);
218  double sinazi = sin(azi);
219  double coslat = cos(lat);
220  double sinlat = sin(lat);
221  double coslon = cos(lon);
222  double sinlon = sin(lon);
223  double uNorth = cosalt * cosazi;
224  double uEast = cosalt * sinazi;
225  double uUp = sinalt;
226  double uRho = - sinlat * uNorth + coslat * uUp;
227  u[0] = coslon * uRho - sinlon * uEast;
228  u[1] = sinlon * uRho + coslon * uEast;
229  u[2] = coslat * uNorth + sinlat * uUp;
230  return;
231 }
232 
233 
234 void XLALComputeDetAMResponseParts(double *armlen, double *xcos, double *ycos, double *fxplus, double *fyplus, double *fxcross, double *fycross, const LALDetector *detector, double ra, double dec, double psi, double gmst)
235 {
236  double X[3]; /* wave frame x axis */
237  double Y[3]; /* wave frame y axis */
238  double Z[3]; /* wave frame z axis (propagation direction) */
239  double U[3]; /* x arm unit vector */
240  double V[3]; /* y arm unit vector */
241  double DU[3][3]; /* single arm response tensor for x arm */
242  double DV[3][3]; /* single arm response tensor for y arm */
243  double gha = gmst - ra; /* greenwich hour angle */
244  double cosgha = cos(gha);
245  double singha = sin(gha);
246  double cosdec = cos(dec);
247  double sindec = sin(dec);
248  double cospsi = cos(psi);
249  double sinpsi = sin(psi);
250  int i, j;
251 
252  /* compute unit vectors specifying the wave frame x, y, and z axes */
253 
254  X[0] = -cospsi * singha - sinpsi * cosgha * sindec;
255  X[1] = -cospsi * cosgha + sinpsi * singha * sindec;
256  X[2] = sinpsi * cosdec;
257  Y[0] = sinpsi * singha - cospsi * cosgha * sindec;
258  Y[1] = sinpsi * cosgha + cospsi * singha * sindec;
259  Y[2] = cospsi * cosdec;
260  Z[0] = -cosgha * cosdec;
261  Z[1] = singha * cosdec;
262  Z[2] = -sindec;
263 
264  switch (detector->type) {
265 
268 
269  /* FIXME: should compute the effect of non-equal arm lengths;
270  * but, for now, just use the mean arm length */
271 
272  *armlen = detector->frDetector.xArmMidpoint
273  + detector->frDetector.yArmMidpoint;
274 
275  /* get the unit vectors along the arms */
276 
278  detector->frDetector.xArmAzimuthRadians,
281 
283  detector->frDetector.yArmAzimuthRadians,
286 
287  /* compute direction cosines for the signal direction relative
288  * to the x-arm and the y-arm */
289 
290  *xcos = *ycos = 0.0;
291  for (i = 0; i < 3; ++i) {
292  *xcos += U[i] * Z[i];
293  *ycos += V[i] * Z[i];
294  }
295 
296  /* compute the single arm response tensors for the x-arm and
297  * y-arm */
298 
299  for (i = 0; i < 3; ++i) {
300  DU[i][i] = 0.5 * U[i] * U[i];
301  DV[i][i] = 0.5 * V[i] * V[i];
302  for (j = i + 1; j < 3; ++j) {
303  DU[i][j] = DU[j][i] = 0.5 * U[i] * U[j];
304  DV[i][j] = DV[j][i] = 0.5 * V[i] * V[j];
305  }
306  }
307 
308  /* compute the beam pattern partial responses for the x-arm and
309  * y-arm */
310 
311  *fxplus = *fxcross = 0.0;
312  *fyplus = *fycross = 0.0;
313  for (i = 0; i < 3; ++i) {
314  double DUX = DU[i][0]*X[0]+DU[i][1]*X[1]+DU[i][2]*X[2];
315  double DUY = DU[i][0]*Y[0]+DU[i][1]*Y[1]+DU[i][2]*Y[2];
316  double DVX = DV[i][0]*X[0]+DV[i][1]*X[1]+DV[i][2]*X[2];
317  double DVY = DV[i][0]*Y[0]+DV[i][1]*Y[1]+DV[i][2]*Y[2];
318  *fxplus += X[i] * DUX - Y[i] * DUY;
319  *fxcross += X[i] * DUY + Y[i] * DUX;
320  *fyplus += X[i] * DVX - Y[i] * DVY;
321  *fycross += X[i] * DVY + Y[i] * DVX;
322  }
323 
324  /* differential interferometer: arm y is subtracted from
325  * arm x */
326  if (detector->type == LALDETECTORTYPE_IFODIFF) {
327  *fyplus *= -1;
328  *fycross *= -1;
329  }
330 
331  break;
332 
334 
335  /* no y-arm */
336 
337  *armlen = 2.0 * detector->frDetector.xArmMidpoint;
338 
340  detector->frDetector.xArmAzimuthRadians,
343 
344  *xcos = *ycos = 0.0;
345  for (i = 0; i < 3; ++i)
346  *xcos += U[i] * Z[i];
347 
348  *fyplus = *fycross = 0.0;
349  XLALComputeDetAMResponse(fxplus, fxcross,
350  (const REAL4(*)[3])detector->response, ra, dec, psi,
351  gmst);
352 
353  break;
354 
356 
357  /* no x-arm */
358 
359  *armlen = 2.0 * detector->frDetector.yArmMidpoint;
360 
362  detector->frDetector.yArmAzimuthRadians,
365 
366  *xcos = *ycos = 0.0;
367  for (i = 0; i < 3; ++i)
368  *ycos += V[i] * Z[i];
369 
370  *fxplus = *fxcross = 0.0;
371  XLALComputeDetAMResponse(fyplus, fycross,
372  (const REAL4(*)[3])detector->response, ra, dec, psi,
373  gmst);
374 
375  break;
376 
377  default:
378 
379  /* FIXME: could handle this situation properly; fur now, just
380  * ignore non long-wavelength-limit effects by setting armlen
381  * to zero; also, pretend that all of the response is
382  * associated with the x-arm */
383 
384  *armlen = *xcos = *ycos = 0.0;
385  *fyplus = *fycross = 0.0;
386  XLALComputeDetAMResponse(fxplus, fxcross,
387  (const REAL4(*)[3])detector->response, ra, dec, psi,
388  gmst);
389 
390  break;
391 
392  }
393 
394  return;
395 }
396 
397 
398 /**
399  * \deprecated Use XLALComputeDetAMResponse() instead.
400  */
401 void LALComputeDetAMResponse(LALStatus * status, LALDetAMResponse * pResponse, const LALDetAndSource * pDetAndSrc, const LIGOTimeGPS * gps)
402 {
403  double fplus, fcross;
404 
405  INITSTATUS(status);
406  ATTATCHSTATUSPTR(status);
407 
408  ASSERT(pResponse != (LALDetAMResponse *) NULL, status, DETRESPONSEH_ENULLOUTPUT, DETRESPONSEH_MSGENULLOUTPUT);
409 
410  ASSERT(pDetAndSrc != NULL, status, DETRESPONSEH_ENULLINPUT, DETRESPONSEH_MSGENULLINPUT);
411 
412  ASSERT(gps != (LIGOTimeGPS *) NULL, status, DETRESPONSEH_ENULLINPUT, DETRESPONSEH_MSGENULLINPUT);
413 
414  /* source coordinates must be in equatorial system */
415  ASSERT(pDetAndSrc->pSource->equatorialCoords.system == COORDINATESYSTEM_EQUATORIAL, status, DETRESPONSEH_ESRCNOTEQUATORIAL, DETRESPONSEH_MSGESRCNOTEQUATORIAL);
416 
418 
419  pResponse->plus = fplus;
420  pResponse->cross = fcross;
421  pResponse->scalar = 0.0; /* not implemented */
422 
423  DETATCHSTATUSPTR(status);
424  RETURN(status);
425 }
426 
427 
428 /**
429  * Computes REAL4TimeSeries containing time series of response amplitudes.
430  * \see XLALComputeDetAMResponse() for more details.
431  */
432 int XLALComputeDetAMResponseSeries(REAL4TimeSeries ** fplus, REAL4TimeSeries ** fcross, const REAL4 D[3][3], const double ra, const double dec, const double psi, const LIGOTimeGPS * start, const double deltaT, const int n)
433 {
434  LIGOTimeGPS t;
435  double gmst;
436  int i;
437  double p, c;
438 
439  *fplus = XLALCreateREAL4TimeSeries("plus", start, 0.0, deltaT, &lalDimensionlessUnit, n);
440  *fcross = XLALCreateREAL4TimeSeries("cross", start, 0.0, deltaT, &lalDimensionlessUnit, n);
441  if(!*fplus || !*fcross) {
444  *fplus = *fcross = NULL;
446  }
447 
448  for(i = 0; i < n; i++) {
449  t = *start;
450  gmst = XLALGreenwichMeanSiderealTime(XLALGPSAdd(&t, i * deltaT));
451  if(XLAL_IS_REAL8_FAIL_NAN(gmst)) {
454  *fplus = *fcross = NULL;
456  }
457  XLALComputeDetAMResponse(&p, &c, D, ra, dec, psi, gmst);
458  (*fplus)->data->data[i] = p;
459  (*fcross)->data->data[i] = c;
460  }
461 
462  return 0;
463 }
464 
465 /**
466  * Computes REAL4TimeSeries containing time series of the full general
467  * metric theory of gravity response amplitudes.
468  * \see XLALComputeDetAMResponseExtraModes() for more details.
469  */
470 int XLALComputeDetAMResponseExtraModesSeries(REAL4TimeSeries ** fplus, REAL4TimeSeries ** fcross, REAL4TimeSeries ** fb, REAL4TimeSeries ** fl, REAL4TimeSeries ** fx, REAL4TimeSeries ** fy, const REAL4 D[3][3], const double ra, const double dec, const double psi, const LIGOTimeGPS * start, const double deltaT, const int n)
471 {
472  LIGOTimeGPS t;
473  double gmst;
474  int i;
475  double p, c, b, l, x, y;
476 
477  *fplus = XLALCreateREAL4TimeSeries("plus", start, 0.0, deltaT, &lalDimensionlessUnit, n);
478  *fcross = XLALCreateREAL4TimeSeries("cross", start, 0.0, deltaT, &lalDimensionlessUnit, n);
479  *fb = XLALCreateREAL4TimeSeries("b", start, 0.0, deltaT, &lalDimensionlessUnit, n);
480  *fl = XLALCreateREAL4TimeSeries("l", start, 0.0, deltaT, &lalDimensionlessUnit, n);
481  *fx = XLALCreateREAL4TimeSeries("x", start, 0.0, deltaT, &lalDimensionlessUnit, n);
482  *fy = XLALCreateREAL4TimeSeries("y", start, 0.0, deltaT, &lalDimensionlessUnit, n);
483 
484  if(!*fplus || !*fcross || !*fb || !*fl || !*fx || !*fy) {
491 
492  *fplus = *fcross = *fb = *fl = *fx = *fy = NULL;
494  }
495 
496  for(i = 0; i < n; i++) {
497  t = *start;
498  gmst = XLALGreenwichMeanSiderealTime(XLALGPSAdd(&t, i * deltaT));
499  if(XLAL_IS_REAL8_FAIL_NAN(gmst)) {
506 
507  *fplus = *fcross = *fb = *fl = *fx = *fy = NULL;
509  }
510  XLALComputeDetAMResponseExtraModes(&p, &c, &b, &l, &x, &y, D, ra, dec, psi, gmst);
511  (*fplus)->data->data[i] = p;
512  (*fcross)->data->data[i] = c;
513  (*fb)->data->data[i] = b;
514  (*fl)->data->data[i] = l;
515  (*fx)->data->data[i] = x;
516  (*fy)->data->data[i] = y;
517  }
518 
519  return 0;
520 }
521 
522 /**
523  * Computes REAL4TimeSeries containing time series of response amplitudes.
524  * \deprecated Use XLALComputeDetAMResponseSeries() instead.
525  */
526 void LALComputeDetAMResponseSeries(LALStatus * status, LALDetAMResponseSeries * pResponseSeries, const LALDetAndSource * pDetAndSource, const LALTimeIntervalAndNSample * pTimeInfo)
527 {
528  /* Want to loop over the time and call LALComputeDetAMResponse() */
529  LALDetAMResponse instResponse;
530  unsigned i;
531 
532  INITSTATUS(status);
533  ATTATCHSTATUSPTR(status);
534 
535  /*
536  * Error-checking assertions
537  */
538  ASSERT(pResponseSeries != (LALDetAMResponseSeries *) NULL, status, DETRESPONSEH_ENULLOUTPUT, DETRESPONSEH_MSGENULLOUTPUT);
539 
540  ASSERT(pDetAndSource != (LALDetAndSource *) NULL, status, DETRESPONSEH_ENULLINPUT, DETRESPONSEH_MSGENULLINPUT);
541 
542  ASSERT(pTimeInfo != (LALTimeIntervalAndNSample *) NULL, status, DETRESPONSEH_ENULLINPUT, DETRESPONSEH_MSGENULLINPUT);
543 
544  /*
545  * Set names
546  */
547  pResponseSeries->pPlus->name[0] = '\0';
548  pResponseSeries->pCross->name[0] = '\0';
549  pResponseSeries->pScalar->name[0] = '\0';
550 
551  strncpy(pResponseSeries->pPlus->name, "plus", LALNameLength);
552  strncpy(pResponseSeries->pCross->name, "cross", LALNameLength);
553  strncpy(pResponseSeries->pScalar->name, "scalar", LALNameLength);
554 
555  /*
556  * Set sampling parameters
557  */
558  pResponseSeries->pPlus->epoch = pTimeInfo->epoch;
559  pResponseSeries->pPlus->deltaT = pTimeInfo->deltaT;
560  pResponseSeries->pPlus->f0 = 0.;
561  pResponseSeries->pPlus->sampleUnits = lalDimensionlessUnit;
562 
563  pResponseSeries->pCross->epoch = pTimeInfo->epoch;
564  pResponseSeries->pCross->deltaT = pTimeInfo->deltaT;
565  pResponseSeries->pCross->f0 = 0.;
566  pResponseSeries->pCross->sampleUnits = lalDimensionlessUnit;
567 
568  pResponseSeries->pScalar->epoch = pTimeInfo->epoch;
569  pResponseSeries->pScalar->deltaT = pTimeInfo->deltaT;
570  pResponseSeries->pScalar->f0 = 0.;
571  pResponseSeries->pScalar->sampleUnits = lalDimensionlessUnit;
572 
573  /*
574  * Ensure enough memory for requested vectors
575  */
576  if(pResponseSeries->pPlus->data->length < pTimeInfo->nSample) {
577  if(lalDebugLevel >= 8)
578  LALInfo(status, "plus sequence too short -- reallocating");
579 
580  TRY(LALSDestroyVector(status->statusPtr, &(pResponseSeries->pPlus->data)), status);
581 
582  TRY(LALSCreateVector(status->statusPtr, &(pResponseSeries->pPlus->data), pTimeInfo->nSample), status);
583 
584  if(lalDebugLevel > 0)
585  printf("pResponseSeries->pPlus->data->length = %d\n", pResponseSeries->pPlus->data->length);
586 
587  }
588 
589  if(pResponseSeries->pCross->data->length < pTimeInfo->nSample) {
590  if(lalDebugLevel >= 8)
591  LALInfo(status, "cross sequence too short -- reallocating");
592 
593  TRY(LALSDestroyVector(status->statusPtr, &(pResponseSeries->pCross->data)), status);
594 
595  TRY(LALSCreateVector(status->statusPtr, &(pResponseSeries->pCross->data), pTimeInfo->nSample), status);
596 
597  }
598 
599  if(pResponseSeries->pScalar->data->length < pTimeInfo->nSample) {
600  if(lalDebugLevel & 0x08)
601  LALInfo(status, "scalar sequence too short -- reallocating");
602 
603  TRY(LALSDestroyVector(status->statusPtr, &(pResponseSeries->pScalar->data)), status);
604 
605  TRY(LALSCreateVector(status->statusPtr, &(pResponseSeries->pScalar->data), pTimeInfo->nSample), status);
606  }
607 
608 
609  /*
610  * Loop to compute each element in time series.
611  */
612 
613  for(i = 0; i < pTimeInfo->nSample; ++i) {
614  LIGOTimeGPS gps = pTimeInfo->epoch;
615  XLALGPSAdd(&gps, i * pTimeInfo->deltaT);
616 
617  TRY(LALComputeDetAMResponse(status->statusPtr, &instResponse, pDetAndSource, &gps), status);
618 
619  pResponseSeries->pPlus->data->data[i] = instResponse.plus;
620  pResponseSeries->pCross->data->data[i] = instResponse.cross;
621  pResponseSeries->pScalar->data->data[i] = instResponse.scalar;
622  }
623 
624  DETATCHSTATUSPTR(status);
625  RETURN(status);
626 }
struct tagLALStatus * statusPtr
Pointer to the next node in the list; NULL if this function is not reporting a subroutine error...
Definition: LALDatatypes.h:954
void XLALComputeDetAMResponseParts(double *armlen, double *xcos, double *ycos, double *fxplus, double *fyplus, double *fxcross, double *fycross, const LALDetector *detector, double ra, double dec, double psi, double gmst)
Definition: DetResponse.c:234
IFO in one-armed mode (X arm)
Definition: LALDetectors.h:229
IFO in one-armed mode (Y arm)
Definition: LALDetectors.h:230
LALSource * pSource
Pointer to LALSource object containing information about the source.
Definition: DetResponse.h:130
SkyPosition equatorialCoords
equatorial coordinates of source, in decimal RADIANS
Definition: DetResponse.h:107
void XLALComputeDetAMResponseExtraModes(double *fplus, double *fcross, double *fb, double *fl, double *fx, double *fy, const REAL4 D[3][3], const double ra, const double dec, const double psi, const double gmst)
An implementation of the detector response for all six tensor, vector and scalar polarisation modes o...
Definition: DetResponse.c:112
#define DETRESPONSEH_ENULLOUTPUT
Output is NULL.
Definition: DetResponse.h:80
REAL4 xArmAzimuthRadians
The angle clockwise from North to the projection of the X arm (or bar&#39;s cylidrical axis) into the lo...
Definition: LALDetectors.h:250
void LALSCreateVector(LALStatus *, REAL4Vector **, UINT4)
REAL4 scalar
Detector response to scalar gravitational radiation (NB: ignored at present – scalar response comput...
Definition: DetResponse.h:143
COMPLEX16 XLALComputeDetArmTransferFunction(double beta, double mu)
Definition: DetResponse.c:202
LALUnit sampleUnits
The physical units of the quantity being sampled.
Definition: LALDatatypes.h:575
Epoch relative to GPS epoch, see LIGOTimeGPS type for more details.
Definition: LALDatatypes.h:458
REAL8 XLALGreenwichMeanSiderealTime(const LIGOTimeGPS *gpstime)
Convenience wrapper, calling XLALGreenwichSiderealTime() with the equation of equinoxes set to 0...
The sky-fixed equatorial coordinate system.
REAL4TimeSeries * pPlus
timeseries of detector response to -polarized gravitational radiation
Definition: DetResponse.h:154
REAL8 latitude
The latitudinal coordinate (in radians), as defined above.
#define XLAL_ERROR(...)
Macro to invoke a failure from a XLAL routine returning an integer.
Definition: XLALError.h:701
Detector structure.
Definition: LALDetectors.h:265
REAL8 vertexLongitudeRadians
The geodetic longitude of the vertex in radians.
Definition: LALDetectors.h:246
int XLALComputeDetAMResponseSeries(REAL4TimeSeries **fplus, REAL4TimeSeries **fcross, const REAL4 D[3][3], const double ra, const double dec, const double psi, const LIGOTimeGPS *start, const double deltaT, const int n)
Computes REAL4TimeSeries containing time series of response amplitudes.
Definition: DetResponse.c:432
This structure aggregates a pointer to a LALDetector and a LALSource.
Definition: DetResponse.h:126
REAL4 xArmAltitudeRadians
The angle up from the local tangent plane of the reference ellipsoid to the X arm (or bar&#39;s cylidric...
Definition: LALDetectors.h:249
int XLALComputeDetAMResponseExtraModesSeries(REAL4TimeSeries **fplus, REAL4TimeSeries **fcross, REAL4TimeSeries **fb, REAL4TimeSeries **fl, REAL4TimeSeries **fx, REAL4TimeSeries **fy, const REAL4 D[3][3], const double ra, const double dec, const double psi, const LIGOTimeGPS *start, const double deltaT, const int n)
Computes REAL4TimeSeries containing time series of the full general metric theory of gravity response...
Definition: DetResponse.c:470
CoordinateSystem system
The coordinate system in which latitude/longitude are expressed.
float REAL4
Single precision real floating-point number (4 bytes).
LIGOTimeGPS * XLALGPSAdd(LIGOTimeGPS *epoch, REAL8 dt)
Adds a double to a GPS time.
Definition: XLALTime.c:131
void XLALComputeDetAMResponse(double *fplus, double *fcross, const REAL4 D[3][3], const double ra, const double dec, const double psi, const double gmst)
An implementation of the detector response formulae in Anderson et al PRD 63 042003 (2001) ...
Definition: DetResponse.c:44
#define LALInfo(statusptr, info)
Definition: LALError.h:119
const LALDetector * pDetector
Pointer to LALDetector object containing information about the detector.
Definition: DetResponse.h:129
REAL4 response[3][3]
The Earth-fixed Cartesian components of the detector&#39;s response tensor .
Definition: LALDetectors.h:268
static double Y(int length, int i)
Maps the length of a window and the offset within the window to the "y" co-ordinate of the LAL docume...
Definition: Window.c:109
IFO in common mode.
Definition: LALDetectors.h:231
void LALComputeDetAMResponse(LALStatus *status, LALDetAMResponse *pResponse, const LALDetAndSource *pDetAndSrc, const LIGOTimeGPS *gps)
Definition: DetResponse.c:401
REAL8 orientation
Orientation angle ( ) of source: counter-clockwise angle -axis makes with a line perpendicular to mer...
Definition: DetResponse.h:108
const LALUnit lalDimensionlessUnit
dimensionless units
Definition: UnitDefs.c:156
void LALComputeDetAMResponseSeries(LALStatus *status, LALDetAMResponseSeries *pResponseSeries, const LALDetAndSource *pDetAndSource, const LALTimeIntervalAndNSample *pTimeInfo)
Computes REAL4TimeSeries containing time series of response amplitudes.
Definition: DetResponse.c:526
REAL8 f0
The heterodyning frequency, in Hertz (zero if not heterodyned).
Definition: LALDatatypes.h:574
REAL8 vertexLatitudeRadians
The geodetic latitude of the vertex in radians.
Definition: LALDetectors.h:247
#define ATTATCHSTATUSPTR(statusptr)
#define RETURN(statusptr)
UINT4 nSample
The total number of samples to be computed.
Definition: DetResponse.h:171
#define INITSTATUS(statusptr)
REAL4 yArmAzimuthRadians
The angle clockwise from North to the projection of the Y arm into the local tangent plane of the re...
Definition: LALDetectors.h:252
LIGOTimeGPS epoch
The start time of the time series.
Definition: LALDatatypes.h:572
LALDetectorType type
The type of the detector (e.g., IFO in differential mode, cylindrical bar, etc.)
Definition: LALDetectors.h:269
void LALSDestroyVector(LALStatus *, REAL4Vector **)
#define DETRESPONSEH_ENULLINPUT
Input is NULL.
Definition: DetResponse.h:79
#define DETATCHSTATUSPTR(statusptr)
#define lalDebugLevel
Definition: LALDebugLevel.h:56
#define XLAL_IS_REAL8_FAIL_NAN(val)
Tests if val is a XLAL REAL8 failure NaN.
Definition: XLALError.h:394
REAL4TimeSeries * XLALCreateREAL4TimeSeries(const CHAR *name, const LIGOTimeGPS *epoch, REAL8 f0, REAL8 deltaT, const LALUnit *sampleUnits, size_t length)
REAL4 * data
Pointer to the data array.
Definition: LALDatatypes.h:150
CHAR name[LALNameLength]
The name of the time series.
Definition: LALDatatypes.h:571
#define DETRESPONSEH_ESRCNOTEQUATORIAL
Source coordinates not in Equatorial system.
Definition: DetResponse.h:81
LAL status structure, see The LALStatus structure for more details.
Definition: LALDatatypes.h:947
REAL4TimeSeries * pCross
timeseries of detector response to -polarized gravitational radiation
Definition: DetResponse.h:155
LIGOTimeGPS epoch
The start time of the time series.
Definition: DetResponse.h:169
This structure encapsulates time and sampling information for computing a LALDetAMResponseSeries.
Definition: DetResponse.h:166
REAL4 xArmMidpoint
The distance to the midpoint of the X arm in meters (unused for bars: set it to zero) ...
Definition: LALDetectors.h:253
This structure encapsulates the detector AM (beam pattern) coefficients for one source at one instanc...
Definition: DetResponse.h:138
REAL4TimeSeries * pScalar
timeseries of detector response to scalar gravitational radiation (NB: not yet implemented.)
Definition: DetResponse.h:156
REAL8 longitude
The longitudinal coordinate (in radians), as defined above.
REAL4 yArmMidpoint
The distance to the midpoint of the Y arm in meters (unused for bars: set it to zero) ...
Definition: LALDetectors.h:254
double complex COMPLEX16
Double-precision floating-point complex number (16 bytes total)
REAL4 yArmAltitudeRadians
The angle up from the local tangent plane of the reference ellipsoid to the Y arm in radians (unused...
Definition: LALDetectors.h:251
REAL4Sequence * data
The sequence of sampled data.
Definition: LALDatatypes.h:576
This structure aggregates together three REAL4TimeSeries objects containing time series of detector A...
Definition: DetResponse.h:151
Time series of REAL4 data, see DATATYPE-TimeSeries types for more details.
Definition: LALDatatypes.h:570
REAL8 deltaT
The time step between samples of the time series in seconds.
Definition: LALDatatypes.h:573
void XLALDestroyREAL4TimeSeries(REAL4TimeSeries *series)
IFO in differential mode.
Definition: LALDetectors.h:228
static void getarm(double u[3], double alt, double azi, double lat, double lon)
Definition: DetResponse.c:213
#define TRY(func, statusptr)
LALFrDetector frDetector
The original LALFrDetector structure from which this was created.
Definition: LALDetectors.h:270
static double psi(double theta, double xi)
Definition: XLALMarcumQ.c:273
REAL4 plus
Detector response to -polarized gravitational radiation.
Definition: DetResponse.h:141
Internal function call failed bit: "or" this with existing error number.
Definition: XLALError.h:463
#define LAL_PI
Archimedes&#39;s constant, pi.
Definition: LALConstants.h:179
#define ASSERT(assertion, statusptr, code, mesg)
UINT4 length
Number of elements in array.
Definition: LALDatatypes.h:149
REAL8 deltaT
The sampling interval , in seconds.
Definition: DetResponse.h:170
REAL4 cross
Detector response to -polarized gravitational radiation.
Definition: DetResponse.h:142