/*--------------main function---------------*/ int main(int argc, char **argv){ const CHAR *fn = __func__; InputParams XLAL_INIT_DECL(inputs); REAL8 srate = 16384.0; /*sample rate defaulted to 16384 */ /* read in command line input args */ ReadInput( &inputs, argc, argv ); LALStatus XLAL_INIT_DECL(status); EphemerisData *edat; if ( (edat = InitEphemeris ( inputs.ephemType, inputs.ephemDir)) == NULL ){ XLALPrintError ( "%s: Failed to init ephemeris data\n", fn ); XLAL_ERROR ( XLAL_EFUNC ); } /*init detector info */ LALDetector *site; if ( ( site = XLALGetSiteInfo ( inputs.det )) == NULL ){ XLALPrintError("%s: Failed to get site-info for detector '%s'\n", fn, inputs.det ); XLAL_ERROR ( XLAL_EFUNC ); } if( inputs.geocentre ){ /* set site to the geocentre */ site->location[0] = 0.0; site->location[1] = 0.0; site->location[2] = 0.0; } struct dirent **pulsars; INT4 n=scandir(inputs.pulsarDir, &pulsars, 0, alphasort); if ( n < 0){ XLALPrintError("scandir failed\n"); XLAL_ERROR(XLAL_EIO); } UINT4 numpulsars = (UINT4)n; UINT4 h=0; CHAR parname[256]; PulsarParameters *pulparams[numpulsars]; for(h=2; h<numpulsars; h++){ if(strstr(pulsars[h]->d_name,".par") == NULL){ free(pulsars[h]); continue; } else{ sprintf(parname,"%s/%s", inputs.pulsarDir, pulsars[h]->d_name); fprintf(stderr, "%s\n", parname); FILE *inject; if (( inject = fopen ( parname, "r" )) == NULL ){ fprintf(stderr,"Error opening file: %s\n", parname); XLAL_ERROR ( XLAL_EIO ); } pulparams[h] = XLALReadTEMPOParFile( parname ); fclose( inject ); } } LIGOTimeGPS epoch; UINT4 ndata; epoch.gpsSeconds = inputs.epoch; epoch.gpsNanoSeconds = 0; ndata = inputs.frDur; REAL8TimeSeries *series=NULL; CHAR out_file[256]; sprintf(out_file, "%s-%s-%d-%d.gwf", inputs.det, inputs.outStr, epoch.gpsSeconds, ndata ); LALFrameH *outFrame = NULL; if ((outFrame = XLALFrameNew( &epoch, (REAL8)ndata, inputs.channel, 1, 0, 0 )) == NULL) { LogPrintf(LOG_CRITICAL, "%s : XLALFrameNew() filed with error = %d.\n", fn, xlalErrno); XLAL_ERROR( XLAL_EFAILED); } if ((series = XLALCreateREAL8TimeSeries( inputs.channel, &epoch, 0., 1./srate,&lalSecondUnit, (int)(ndata*srate) )) == NULL) { XLAL_ERROR( XLAL_EFUNC ); } UINT4 counter=0; for (counter = 0; counter < series->data->length; counter++) series->data->data[counter] = 0; /*** Read Pulsar Data ***/ for (h=0; h < numpulsars; h++){ if(strstr(pulsars[h]->d_name,".par")==NULL){ free(pulsars[h]); continue; } else{ PulsarSignalParams XLAL_INIT_DECL(params); /* set signal generation barycenter delay look-up table step size */ params.dtDelayBy2 = 10.; /* generate table every 10 seconds */ if (( params.pulsar.spindown = XLALCreateREAL8Vector(1)) == NULL ){ XLALPrintError("Out of memory"); XLAL_ERROR ( XLAL_EFUNC ); } INT4 dtpos = 0; if ( PulsarCheckParam(pulparams[h], "POSEPOCH") ) dtpos = epoch.gpsSeconds - (INT4)PulsarGetREAL8Param(pulparams[h], "POSEPOCH"); else dtpos = epoch.gpsSeconds - (INT4)PulsarGetREAL8Param(pulparams[h], "PEPOCH"); REAL8 ra = 0., dec = 0.; if ( PulsarCheckParam( pulparams[h], "RAJ" ) ) { ra = PulsarGetREAL8Param( pulparams[h], "RAJ" ); } else if ( PulsarCheckParam( pulparams[h], "RA" ) ){ ra = PulsarGetREAL8Param( pulparams[h], "RA" ); } else{ XLALPrintError("No right ascension found"); XLAL_ERROR ( XLAL_EFUNC ); } if ( PulsarCheckParam( pulparams[h], "DECJ" ) ) { dec = PulsarGetREAL8Param( pulparams[h], "DECJ" ); } else if ( PulsarCheckParam( pulparams[h], "DEC" ) ){ dec = PulsarGetREAL8Param( pulparams[h], "DEC" ); } else{ XLALPrintError("No declination found"); XLAL_ERROR ( XLAL_EFUNC ); } params.pulsar.position.latitude = dec + (REAL8)dtpos * PulsarGetREAL8ParamOrZero(pulparams[h], "PMDEC"); params.pulsar.position.longitude = ra + (REAL8)dtpos * PulsarGetREAL8ParamOrZero(pulparams[h], "PMRA") / cos(params.pulsar.position.latitude); params.pulsar.position.system = COORDINATESYSTEM_EQUATORIAL; REAL8Vector *fs = PulsarGetREAL8VectorParam(pulparams[h], "F"); if ( fs->length == 0 ){ XLALPrintError("No frequencies found"); XLAL_ERROR ( XLAL_EFUNC ); } params.pulsar.f0 = 2.*fs->data[0]; if ( fs->length > 1 ){ params.pulsar.spindown->data[0] = 2.*fs->data[1]; } if (( XLALGPSSetREAL8(&(params.pulsar.refTime), PulsarGetREAL8Param(pulparams[h], "PEPOCH")) ) == NULL ) XLAL_ERROR ( XLAL_EFUNC ); params.pulsar.psi = PulsarGetREAL8ParamOrZero(pulparams[h], "PSI"); params.pulsar.phi0 = PulsarGetREAL8ParamOrZero(pulparams[h], "PHI0"); REAL8 cosiota = PulsarGetREAL8ParamOrZero(pulparams[h], "COSIOTA"); REAL8 h0 = PulsarGetREAL8ParamOrZero(pulparams[h], "H0"); params.pulsar.aPlus = 0.5 * h0 * (1. + cosiota * cosiota ); params.pulsar.aCross = h0 * cosiota; /*Add binary later if needed!*/ params.site = site; params.ephemerides = edat; params.startTimeGPS = epoch; params.duration = ndata; params.samplingRate = srate; params.fHeterodyne = 0.; REAL4TimeSeries *TSeries = NULL; LALGeneratePulsarSignal( &status, &TSeries, ¶ms ); if (status.statusCode){ fprintf(stderr, "LAL Routine failed!\n"); XLAL_ERROR (XLAL_EFAILED); } UINT4 i; for (i=0; i < TSeries->data->length; i++) series->data->data[i] += TSeries->data->data[i]; XLALDestroyREAL4TimeSeries(TSeries); XLALDestroyREAL8Vector(params.pulsar.spindown); } } if (XLALFrameAddREAL8TimeSeriesProcData(outFrame,series)){ LogPrintf(LOG_CRITICAL, "%s : XLALFrameAddREAL8TimeSeries() failed with error = %d.\n",fn,xlalErrno); XLAL_ERROR(XLAL_EFAILED); } CHAR OUTFILE[256]; sprintf(OUTFILE, "%s/%s", inputs.outDir, out_file); if ( XLALFrameWrite(outFrame, OUTFILE)){ LogPrintf(LOG_CRITICAL, "%s : XLALFrameWrite() failed with error = %d.\n", fn, xlalErrno); XLAL_ERROR( XLAL_EFAILED ); } XLALFrameFree(outFrame); XLALDestroyREAL8TimeSeries( series ); return 0; }
int main( void ){ INT4 i = 0; PulsarParameters *pars = NULL; /* output par file */ FILE *fp = NULL; if ( (fp = fopen(PARFILE, "w")) == NULL ){ XLAL_ERROR( XLAL_EIO, "Error... problem writing parfile!\n" ); } for( i=0; i<NUMPARS; i++ ){ if ( p[i].sigma != NULL ){ /* write out value and error */ fprintf(fp, "%s\t%s\t%s\t%s\n", p[i].name, p[i].val, p[i].fitFlag, p[i].sigma); } else{ /* just write out value */ fprintf(fp, "%s\t%s\n", p[i].name, p[i].val); } } fclose(fp); /* read in par file */ pars = XLALReadTEMPOParFileNew( PARFILE ); /* check read-in parameters against originals */ for ( i=0; i<NUMPARS; i++ ){ CHAR outval[256]; CHAR outsigma[256]; CHAR outfitflag[256]; /* see if the parameter can be split into an alphabet and numerical value (as is the case * for e.g. FB, which is held as a vector */ CHAR *namecopy = NULL, *token = NULL; namecopy = XLALStringDuplicate( p[i].name ); token = strtok(namecopy, delimiters); /* get part of the name before the numerical delimiter */ if ( !PulsarCheckParam( pars, p[i].name ) && !PulsarCheckParam( pars, token ) ){ XLAL_ERROR( XLAL_EFAILED, "Error... parameter %s does not exist in the read-in data!\n", p[i].name ); } /* value is a vector */ if ( PulsarCheckParam( pars, token ) && PulsarGetParamType( pars, token ) == PULSARTYPE_REAL8Vector_t ){ /* get value and convert to string */ if ( !strchr( p[i].valcheck, 'e' ) ){ /* number doesn't contain an exponent */ sprintf(outval, "%.5lf", PulsarGetREAL8VectorParamIndividual(pars, p[i].name)); } else{ /* number does contain an exponent */ sprintf(outval, "%.5le", PulsarGetREAL8VectorParamIndividual(pars, p[i].name)); } if ( p[i].sigmacheck != NULL ){ /* get error and convert to string */ if ( !strchr( p[i].sigmacheck, 'e' ) ){ /* number doesn't contain an exponent */ sprintf(outsigma, "%.5lf", PulsarGetREAL8VectorParamErrIndividual(pars, p[i].name)); } else{ sprintf(outsigma, "%.5le", PulsarGetREAL8VectorParamErrIndividual(pars, p[i].name)); } UINT4 *fitFlag = PulsarGetParamFitFlag(pars, token); if ( fitFlag[atoi(p[i].name+strlen(token))] == 1 ){ sprintf(outfitflag, "1"); } else if ( fitFlag[atoi(p[i].name+strlen(token))] == 0 ){ sprintf(outfitflag, " "); } else { XLAL_ERROR( XLAL_EFAILED, "Error... fit flag incorrect for %s.\n", p[i].name ); } } } /* value is a double */ else if( PulsarGetParamType( pars, p[i].name ) == PULSARTYPE_REAL8_t ){ /* get value and convert to string */ if ( !strchr( p[i].valcheck, 'e' ) ){ /* number doesn't contain an exponent */ sprintf(outval, "%.5lf", PulsarGetREAL8Param(pars, p[i].name)); } else{ /* number does contain an exponent */ sprintf(outval, "%.5le", PulsarGetREAL8Param(pars, p[i].name)); } if ( p[i].sigmacheck != NULL ){ /* get error and convert to string */ if ( !strchr( p[i].sigmacheck, 'e' ) ){ /* number doesn't contain an exponent */ sprintf(outsigma, "%.5lf", PulsarGetREAL8ParamErr(pars, p[i].name)); } else{ sprintf(outsigma, "%.5le", PulsarGetREAL8ParamErr(pars, p[i].name)); } UINT4 *fitFlag = PulsarGetParamFitFlag(pars, p[i].name); if ( fitFlag[0] == 1 ){ sprintf(outfitflag, "1"); } else if ( fitFlag[0] == 0 ){ sprintf(outfitflag, " "); } else { XLAL_ERROR( XLAL_EFAILED, "Error... fit flag incorrect for %s.\n", p[i].name ); } } } /* value is a string */ else if ( PulsarGetParamType( pars, p[i].name ) == PULSARTYPE_string_t ) { CHAR *out = XLALStringDuplicate(PulsarGetStringParam(pars, p[i].name)); sprintf(outval, "%s", out); XLALFree(out); } /* compare returned value with input value */ if ( strcmp(outval, p[i].valcheck) != 0 ){ XLAL_ERROR( XLAL_EFAILED, "Error... parameter %s does not match input (%s cf. %s)!\n", p[i].name, outval, p[i].valcheck ); } if( p[i].sigma != NULL ){ /* compare returned value with input value */ if ( strcmp(outsigma, p[i].sigmacheck) != 0 ){ XLAL_ERROR( XLAL_EFAILED, "Error... parameter sigma %s does not match input (%s cf. %s)!\n", p[i].name, outsigma, p[i].sigmacheck ); } if ( strcmp(outfitflag, p[i].fitFlag) != 0 ){ XLAL_ERROR( XLAL_EFAILED, "Error... parameter %s fit flag does not match input!\n", p[i].fitFlag ); } } XLALFree( namecopy ); } /* remove par file */ if ( remove(PARFILE) != 0 ){ XLAL_ERROR( XLAL_EIO, "Error... problem removing parfile!\n" ); } fprintf(stderr, "Successfully read in .par file values!\n"); PulsarFreeParams( pars ); LALCheckMemoryLeaks(); return XLAL_SUCCESS; }