/***************************************************************** * Name of source file: initTslice2.c * Version of source file: 2.3.7 * Parent component: j_ima_iros * Programmer: Niels J. Westergaard * Affiliation: National Space Institute/DTU * Contact address: Juliane Maries Vej 30 * DK-2100 Copenhagen, Denmark * Phone: +45 35325705 FAX: +45 35362475 * njw@space.dtu.dk * Purpose: Create a shadowgram from the event list * Origin date: 071029 * Update history: 2.3.7 100308 First version * Type of element: Function * Part of package: j_ima_iros * Name: initTslice2 *****************************************************************/ #include int initTslice2( timeslices *Tslice, Auxdata auxdata, Events *bufevents, struct log_con *logger, int chatter, int status) { int flag = 0, i5 = 5; long i = 0; long i1 = 0; long i2 = 0; long idx = 0; float gain_sum = 0.0; float deadtime_sum = 0.0; float greyfilter_sum = 0.0; static float gain1=22.0, gain2=22.0; double ontime = 0.0; double time_step = 0.0; double time = 0.0; double *dblBuffer = NULL; double IJDstart; double IJDstop; static int func_local = -1, local_calling_num=0; clock_t calling_TT, TT_now; if( status != ISDC_OK ) return status; /* function entry timing start */ if (func_local < 0) { func_local = logger->func_used_num++; logger->func_calling_num = 0; sprintf(logger->func_names[func_local], "initTslice2"); if (logger->func_used_num > J_IIR_MAX_NUM_FUNCS) { sprintf(logger->tmptxt, "func_used_num > %d", J_IIR_MAX_NUM_FUNCS); logger->logstat = logprint(logger, logger->tmptxt, 6); logger->func_used_num = J_IIR_MAX_NUM_FUNCS; } } logger->func_calls[func_local]++; logger->TT = TT_now = clock(); calling_TT = logger->TT; local_calling_num = logger->func_calling_num; logger->func_calling_num = func_local; /* function entry timing end */ i5 = logger->i5; /* * Initialize Tslice and then fill the gain, deadtime and DontUseFlag * assuming that the DAL3JEMX event selection has taken the GTIs * properly into account * Tslice->Tstart and ->Tstop has been initialized in 'main' */ IJDstart = Tslice->Tstart; IJDstop = Tslice->Tstop; Tslice->onTime = 0.0; Tslice->avGain = 22.0; Tslice->gainclass = 0; Tslice->deadc = 0.0; Tslice->exposure = 0.0; Tslice->avCorf = 0.0; Tslice->srcWCounts = 0.0; Tslice->srcCounts = 0.0; Tslice->bkgWCounts = 0.0; Tslice->bkgCounts = 0.0; Tslice->Events = 0.0; Tslice->wEvents = 0.0; Tslice->flux = 0.0; Tslice->fErr = 0.0; Tslice->burstflag = 0; Tslice->dontUse = 0; /* Find event indices for first and last event in the defined interval */ status = j_iir_bin_search( nEvents, dblBuffer, IJDstart, &i1, status ); status = j_iir_bin_search( nEvents, dblBuffer, IJDstop, &i2, status ); /* * Loop between the indices and update the shadowgram and calculate the average gain */ gain_sum = 0.0; if( i2 >= i1 && bufevents[i1].ijd >= IJDstart && bufevents[i2].ijd <= IJDstop ) { /* fine! */ shd_adhoc.nEvs = 0; for( i = i1; i <= i2; i++ ) { if( bufevents[i].pi < chanMin || bufevents[i].pi > chanMax ) continue; shd_adhoc.nEvs++; idx = bufevents[i].rawx + bufevents[i].rawy * J_IIR_DIM_STDSHDG; shd_adhoc.shd[idx]++; gain_sum += bufevents[i].gain; } shd_adhoc.gain = gain_sum / shd_adhoc.nEvs; } else { shd_adhoc.nEvs = 0; shd_adhoc.gain = 22.9999; /* presumably not used */ } /* * Calculate the ontime in seconds for the time-interval */ if( auxdata.n_GTI > 0 ) { /* skip this step if no GTIs were found */ /* step through the GTIs */ ontime = 0.0; for( i = 0; i < auxdata.n_GTI; i++ ) { /* skip this interval if GTI_IJDstop comes before IJDstart */ if( auxdata.GTI_IJDstop[i] < IJDstart ) continue; /* skip this interval and the following if GTI_IJDstart comes after IJDstop */ if( auxdata.GTI_IJDstart[i] > IJDstop ) break; /* The GTI interval starts before IJDstart */ if( auxdata.GTI_IJDstart[i] < IJDstart ) { /* and the GTI interval stops before IJDstop */ if( auxdata.GTI_IJDstop[i] < IJDstop ) { ontime = (auxdata.GTI_IJDstop[i] - IJDstart)/J_IIR_IJDSECOND; } else { /* the GTI interval covers the entire period */ ontime = (IJDstop - IJDstart)/J_IIR_IJDSECOND; break; } } else { /* The GTI interval start after IJDstart */ /* and the GTI interval stops before IJDstop */ if( auxdata.GTI_IJDstop[i] < IJDstop ) { ontime += (auxdata.GTI_IJDstop[i] - auxdata.GTI_IJDstart[i])/J_IIR_IJDSECOND; } else { /* the GTI interval stops after IJDstop */ ontime += (IJDstop - auxdata.GTI_IJDstart[i])/J_IIR_IJDSECOND; break; } } } } else { /* disregard the GTI information since it is not there */ ontime = (IJDstop - IJDstart)/J_IIR_IJDSECOND; } if (ontime == 0.0) {Tslice->dontUse = 1; goto exittrace;} Tslice->onTime = ontime; Tslice->exposure = ontime; /* preliminary value */ if( ontime > 0.0 ) { /* Start further analysis */ /* * Get the time averaged deadtime (or rather dead-fraction) */ /* Set the timestep to the smaller of 1 s and time-interval/10 */ /* (It is silently assumed that GTIs are of so long durations that we don't need to consider those for the time step.) */ time_step = IJDstop - IJDstart; if( time_step > 10.*J_IIR_IJDSECOND ) { time_step = J_IIR_IJDSECOND; } else if( time_step > J_IIR_IJDSECOND ) { time_step = 0.1*J_IIR_IJDSECOND; } else { time_step *= 0.1; } sprintf(logger->tmptxt, "BCHK: %1d %3d rms: %f, flux: %f, bflag: %1d, signf: %f", PIFmap->BURST_FLAG, ntime, rms[1], Tslice[ntime].wflux[allsrc], Tslice[ntime].burstflag, signif); logger->logstat = logprint(logger, logger->tmptxt, i5); RILlogMessage( NULL, Log_0, "##5111## time_step (in s) : %8.5lf", time_step/J_IIR_IJDSECOND); RILlogMessage( NULL, Log_0, "##5112## auxdata.n_deadTimes = %d", auxdata.n_deadTimes ); RILlogMessage( NULL, Log_0, "##5113## auxdata.n_greyValues = %d", auxdata.n_greyValues ); deadtime_sum = 0.0; greyfilter_sum = 0.0; i2 = 0; for( time = IJDstart; time <= IJDstop; time += time_step ) { /* are we in a GTI ? - else skip */ flag = 0; if( auxdata.n_GTI > 0 ) { for( i = 0; i < auxdata.n_GTI; i++ ) { if( auxdata.GTI_IJDstart[i] <= time && auxdata.GTI_IJDstop[i] >= time ) { RILlogMessage( NULL, Log_0,"##5114## time = %14.8lf inside GTI #%ld", time, i ); flag = 1; break; } } } else { RILlogMessage( NULL, Log_0,"##5115## Disregarding GTI info" ); flag = 1; /* disregard the GTI information */ } if( !flag ) continue; i2++; /* count the number of accepted times */ /* what dead time applies ? */ status = j_iir_bin_search( (long) auxdata.n_deadTimes, auxdata.deadIJD, time, &i, status ); if( chatter > J_CHATTY_VERBOSE ) RILlogMessage( NULL, Log_0, "##5116## auxdata.deadIJD[%ld] = %14.8lf, dead = %7.4lf", i, auxdata.deadIJD[i], auxdata.deadTime[i] ); deadtime_sum += auxdata.deadTime[i]; /* what greyfilter applies ? */ status = j_iir_bin_search( auxdata.n_greyValues, auxdata.greyTimes, time, &i, status ); if( chatter > J_CHATTY_VERBOSE ) RILlogMessage( NULL, Log_0, "##5117## auxdata.greyTimes[%ld] = %14.8lf, grey = %d", i, auxdata.greyTimes[i], auxdata.greyValues[i] ); greyfilter_sum += (auxdata.greyValues[i] + 1)/32.; } shd_adhoc.exposure *= (1.-deadtime_sum/i2)*(greyfilter_sum/i2); } RILlogMessage( NULL, Log_0, "##5118## Final exposure : %9.3lf", shd_adhoc.exposure); exittrace: if (logger->trace) traces(func_local, 99, logger); /* function exit timing start */ TT_now = clock(); logger->func_times[func_local] += difftime(TT_now, logger->TT); logger->TT = TT_now; logger->func_calling_num = local_calling_num; /* function exit timing end */ return status; } /***************************************************************** * Name of source file: initTslice.c * Version of source file: 1.0.0 * Parent component: None * Programmer: Niels Lund * Affiliation: Danish National Space Center * Contact address: Juliane Maries Vej 30 * DK-2100 Copenhagen, Denmark * Phone: +45 35325716 FAX: +45 35362475 * nl@dnsc.dk * Purpose: initialize one member of the Tslices structure * Origin date: 20100226 * Update history: 1.0.0 20100226 First version * Type of element: Function * Part of package: j_ima_iros * Name: initTslice * Required components: DAL, PIL, RIL, DAL3GEN *****************************************************************/ #include /* QQQQQQ */ int initTslice( timeslices Tslice, struct log_con *logger, int status ) { int i, i5; static int func_local = -1, local_calling_num=0; clock_t calling_TT, TT_now; /* function entry timing start */ if (func_local < 0) { func_local = logger->func_used_num++; sprintf(logger->func_names[func_local], "initTslice"); if (logger->func_used_num > 199) { sprintf(logger->tmptxt, "func_used_num > 199"); logger->logstat = logprint(logger, logger->tmptxt, 6); logger->func_used_num = 199; } } logger->func_calls[func_local]++; logger->TT = TT_now = clock(); calling_TT = logger->TT; local_calling_num = logger->func_calling_num; logger->func_calling_num = func_local; /* function entry timing end */ if (logger->trace) traces(func_local, 0, logger); i5 = logger->i5; Tslice.Tstart = 0.0; Tslice.Tstop = 0.0; Tslice.onTime = 0.0; Tslice.avGain = 22.0; Tslice.gainclass = 0; Tslice.deadc = 0.0; Tslice.exposure = 0.0; Tslice.avCorf = 0.0; Tslice.srcWCounts = 0.0; Tslice.srcCounts = 0.0; Tslice.bkgWCounts = 0.0; Tslice.bkgCounts = 0.0; Tslice.Events = 0.0; Tslice.wEvents = 0.0; Tslice.flux = 0.0; Tslice.fErr = 0.0; Tslice.burstflag = 0; /* typedef struct { stucture containing information about one time bin (reused for every user-energy) double Tstart; ijd-scale (kept constant for all energy bins) double Tstop; ijd-scale (kept constant for all energy bins) float onTime; ijd-units (kept constant for all energy bins) float avGain; gain-units (kept constant for all energy bins) int gainclass; dimensionless (kept constant for all energy bins) float deadc; dimensionless (may vary between energy bins) float exposure; ijd-units (may vary between energy bins) float avCorf; dimensionless (may vary between energy bins) float srcWCounts; weighted sum of counts for high-illum pixels float srcCounts; sum of counts for high-illum pixels float bkgWCounts; weighted sum of counts for low-illum pixels float bkgCounts; sum of counts for low-illum pixels float Events; all valid events in time bin float wEvents; all valid events passing E-selection in time bin float flux; final calculated flux - at current energy float fErr; final calculated flux-error - at current energy } timeslices; */ slut: ; if (logger->trace) traces(func_local, 99, logger); /* function exit timing start */ TT_now = clock(); logger->func_times[func_local] += difftime(TT_now, logger->TT); logger->TT = TT_now; logger->func_calling_num = local_calling_num; /* function exit timing end */ return( status ); exittrace: /* exittrace */ if (logger->trace) traces(func_local, 1000-status, logger); goto slut; }