3 #define EXIT(ERR) gSystem->Exit(ERR) // better exit handling for ROOT stuff
21 {cout <<
"cwb1G::Init - Error : analysis must be 1G" << endl;
EXIT(1);}
23 {cout <<
"cwb1G::Init - Error : if analysis=1G then stage must be CWB_STAGE_FULL" << endl;
EXIT(1);}
29 double dt_max = 1./rate_min;
31 cout <<
"cwb1G::Init - Error : lagStep is not a multple of max-time-resolution" << endl;
32 cout <<
"lagStep(sec) : " <<
cfg.
lagStep <<
"\t max dt(sec) : " << dt_max << endl << endl;
57 std::vector<wavearray<double> >
y;
68 {cout <<
"cwb1G::ReadData - Error opening root file : " <<
jname << endl;
EXIT(1);}
69 TDirectory* cdstrain = (TDirectory*)
jfile->Get(
"strain");
70 if(cdstrain==NULL) cdstrain =
jfile->mkdir(
"strain");
78 {cout <<
"cwb1G::ReadData - input rate from frame " << x.
rate()
79 <<
" do not match the one defined in config : " <<
cfg.
inRate << endl;
EXIT(1);}
98 cdstrain->cd();
pTF[
i]->Write(
ifo[i]);
103 if(i>0 && xstart != x.
start()) {
104 cout <<
"cwb1G::ReadData - Error : ifo noise data not synchronized" << endl;
105 cout <<
ifo[
i] <<
" " << x.
start() <<
" != " <<
ifo[0] <<
" " << xstart << endl;
108 if(i>0 && xrate != x.
rate()) {
109 cout <<
"cwb1G::ReadData - Error : ifo noise data have different rates" << endl;
110 cout <<
ifo[
i] <<
" " << x.
rate() <<
" != " <<
ifo[0] <<
" " << xrate << endl;
115 TDirectory* cdmdc = (TDirectory*)
jfile->Get(
"mdc");
116 if(cdmdc==NULL) cdmdc =
jfile->mkdir(
"mdc");
123 {cout <<
"cwb1G::ReadData - input rate from frame " << x.
rate()
124 <<
" do not match the one defined in config : " <<
cfg.
inRate << endl;
EXIT(1);}
149 cdmdc->cd();wM.Write(
ifo[i]);
153 if(xstart != x.
start()) {
154 cout <<
"cwb1G::ReadData - Error : mdc/noise data with different start time" << endl;
155 printf(
"start time : noise = %10.6f - mdc = %10.6f\n",xstart,x.
start());
158 if(xrate != x.
rate()) {
159 cout <<
"cwb1G::ReadData - Error : mdc/noise data with different rate" << endl;
160 printf(
"rate : noise = %10.6f - mdc = %10.6f\n",xrate,x.
rate());
163 if(xsize != x.
size()) {
164 cout <<
"cwb1G::ReadData - Error : mdc/noise data with different buffer size" << endl;
165 printf(
"buffer size : noise = %d - mdc = %lu\n",xsize,x.
size());
174 TDirectory* cdmdc = (TDirectory*)
jfile->Get(
"mdc");
175 if(cdmdc==NULL) cdmdc =
jfile->mkdir(
"mdc");
178 std::vector<double> mdcFactor;
179 for (
int k=0;
k<(
int)
pD[0]->ISNR.size();
k++) {
182 for(
int i=0;
i<
nIFO;
i++) snr+=
pD[0]->ISNR.data[
k];
184 for(
int i=0;
i<
nIFO;
i++) snr+=
pD[
i]->ISNR.data[
k];
187 if(snr>0) mdcFactor.push_back(1./snr);
else mdcFactor.push_back(0.);
190 size_t K = mdcFactor.size();
191 for(
int k=0;
k<
K;
k++) {
192 if(mdcFactor[
k]) cout << k <<
" mdcFactor : " << mdcFactor[
k] << endl;
203 cdmdc->cd();wM.Write(
ifo[i]);
210 for(
int k=0;
k<(
int)K;
k++) {
211 int ilog[5] = {1,3,12,13,14};
212 for(
int l=0;
l<5;
l++) {
213 double mfactor =
l<2 ? mdcFactor[
k] : mdcFactor[
k]*mdcFactor[
k];
241 TDirectory* cdcstrain=NULL;
251 {cout <<
"cwb1G::DataConditioning - Error : file " <<
jname <<
" not found" << endl;
EXIT(1);}
261 int nshift =
int(factor*pWS->
rate());
265 int jstart = nshift<0 ? -nshift : 0;
266 int jstop = nshift<0 ? pWS->
size() : pWS->
size()-nshift;
267 for(
int j=jstart;
j<jstop;
j++) wM.
data[
j+nshift] = pWS->
data[
j];
271 double tshift=nshift/pWS->
rate();
276 for(
int k=0;
k<(
int)
pD[i]->IWFP.size();
k++) {
281 for(
int k=0;
k<(
int)
pD[i]->TIME.size();
k++)
pD[i]->TIME[
k]+=tshift;
285 vector<TString>
ifos(nIFO);
328 cout<<
"After "<<
ifo[
i]<<
" data conditioning"<<endl;
353 vector<TString> delObjList;
378 cout<<
"live time in zero lag: "<<TL<<endl<<endl;
383 {cout <<
"cwb1G::Coherence - Error : file " <<
jname <<
" not found" << endl;
EXIT(1);}
386 char sfactor[8];
sprintf(sfactor,
"%d",ifactor);
403 cout<<
"pixel threshold in units of noise rms: "<<Ao<<endl;
404 cout<<
"2 OR threshold in units of noise var: "<<
cfg.
x2or*Ao*Ao<<endl;
418 wc.
write(
jfile,
"coherence",
"clusters",0,cycle);
419 wc.
write(
jfile,
"coherence",
"clusters",-1,cycle);
420 cout<<wc.
csize()<<
"|"<<wc.
size()<<
" ";cout.flush();
432 char sfactor[8];
sprintf(sfactor,
"%d",ifactor);
454 {cout <<
"cwb1G::SuperCluster - Error : file " <<
jname <<
" not found" << endl;
EXIT(1);}
462 wc.
read(
jfile,
"coherence",
"clusters",0,cycle);
470 cout<<m<<
"|"<<pwc->
size()<<
" ";
480 vector<TString> delObjList;
515 cout<<
"events in the buffer: "<<
NET.
events()<<
"\n";
521 cout<<
"dump ced into "<<
jname<<
"\n";
522 jfile =
new TFile(jname,
"UPDATE");
524 {cout <<
"cwb1G::Likelihood - Error : file " << jname <<
" not found" << endl;
EXIT(1);}
525 TDirectory* cdced = NULL;
526 cdced = (TDirectory*)
jfile->Get(
"ced");
527 if(cdced == NULL) cdced =
jfile->mkdir(
"ced");
530 cout<<
"dump ced into "<<ced_dir<<
"\n";
544 {cout <<
"cwb1G::Likelihood - Error : file " <<
jname <<
" not found" << endl;
EXIT(1);}
char channelNamesMDC[NIFO_MAX][50]
CWB_JOBF_OPTIONS jobfOptions
CWB::frame fr[2 *NIFO_MAX]
virtual size_t size() const
size_t write(const char *file, int app=0)
double dTau
maximum time delay
void white(double dT=0, int wtype=1, double offset=0., double stride=0.)
what it does: see algorithm description in wseries.hh
static size_t GetProcInfo(bool mvirtual=true)
bool Likelihood(int ifactor, char *ced_dir, netevent *output=NULL, TTree *net_tree=NULL, char *outDump=NULL)
void FileGarbageCollector(TString ifName, TString ofName="", vector< TString > delObjList=vector< TString >())
std::vector< netcluster > wc_List
printf("total live time: non-zero lags = %10.1f \n", liveTot)
bool singleDetector
used for the stage stuff
void set2or(double p)
param: threshold
WSeries< double > * pTF[NIFO_MAX]
pointers to detectors
std::vector< double > * getmdcTime()
size_t setsnr(wavearray< double > &, std::vector< double > *, std::vector< double > *, double=5., double=8.)
fprintf(stdout,"start=%f duration=%f rate=%f\n", x.start(), x.size()/x.rate(), x.rate())
virtual void rate(double r)
virtual WSeries< float > variability(double=0., double=-1., double=-1.)
param: first - time window to calculate normalization constants second - low frequency boundary for c...
TString iname
stage benchmark
detector * pD[NIFO_MAX]
noise variability
size_t setIndexMode(size_t=0)
TFile * jfile
output root file
cout<< endl;cout<< "ts size = "<< ts.size()<< " ts rate = "<< ts.rate()<< endl;tf.Forward(ts, wdm);int levels=tf.getLevel();cout<< "tf size = "<< tf.size()<< endl;double dF=tf.resolution();double dT=1./(2 *dF);cout<< "rate(hz) : "<< RATE<< "\t layers : "<< nLAYERS<< "\t dF(hz) : "<< dF<< "\t dT(ms) : "<< dT *1000.<< endl;int itime=TIME_PIXEL_INDEX;int ifreq=FREQ_PIXEL_INDEX;int index=(levels+1)*itime+ifreq;double time=itime *dT;double freq=(ifreq >0)?ifreq *dF:dF/4;cout<< endl;cout<< "PIXEL TIME = "<< time<< " sec "<< endl;cout<< "PIXEL FREQ = "<< freq<< " Hz "<< endl;cout<< endl;wavearray< double > x
std::vector< wavearray< double > * > RWFP
void add(const wavearray< DataType_t > &, int=0, int=0, int=0)
virtual double ReadData(double mdcShift, int ifactor)
void DataConditioning(int ifactor)
void SetOptions(int simulation, double rho, double inRate, bool useSparse=false, char *gtype=const_cast< char * >("png"), int paletteId=0)
virtual size_t supercluster(char atype, double S, bool core)
param: statistic: E - excess power, L - likelihood param: selection threshold T for likelihood cluste...
virtual void start(double s)
std::vector< double > mdcTime
long coherence(double, double=0., double=0.)
param: threshold on lognormal pixel energy (in units of noise rms) param: threshold on total pixel en...
network NET
pointers to WSeries
size_t cluster(int kt, int kf)
param: time gap in pixels return: number of reconstructed clusters
void Forward(size_t k)
param: number of steps
void Coherence(int ifactor)
double setVeto(double=5.)
param: time window around injections
void bandPass1G(double f1=0., double f2=0.)
void setDelayIndex(double rate)
param: MDC log file
std::vector< std::string > mdcList
int getLayer(wavearray< DataType_t > &w, double n)
param: n - layer number
double dataShift[NIFO_MAX]
int Write(double factor, bool fullCED=true)
char channelNamesRaw[NIFO_MAX][50]
friend void CWB_Plugin(TFile *jfile, CWB::config *, network *, WSeries< double > *, TString, int)
COHERENCE.
void SetChannelName(char *chName)
void PrintStageInfo(CWB_STAGE stage, TString comment, bool out=true, bool log=true, TString fname="")
char jname[1024]
job file object
std::vector< wavearray< double > * > IWFP
size_t cpf(const netcluster &, bool=false, int=0)
long likelihood(char='E', double=sqrt(2.), int=0, size_t=0, int=-1, bool=false)
#define EXPORT(TYPE, VAR, CMD)
WSeries< double > * getTFmap()
param: no parameters
void readFrames(char *filename, char *channel, wavearray< double > &w)
virtual void lprFilter(double, int=0, double=0., double=0.)
void Forward(int n=-1)
param: wavelet - n is number of steps (-1 means full decomposition)
Meyer< double > S(1024, 2)
netcluster * getwc(size_t n)
param: delay index
sprintf(tfres,"(1/%g)x(%g) (sec)x(Hz)", 2 *df, df)
void setDelayFilters(detector *=NULL)
double factors[FACTORS_MAX]
size_t netcut(double, char='L', size_t=0, int=1)
param: threshold param: minimum cluster size processed by the corrcut param: cluster type return: num...
size_t read(const char *)
virtual void resize(unsigned int)
void Inverse(int n=-1)
param: n - number of steps (-1 means full reconstruction)
void putLayer(wavearray< DataType_t > &, double n)
param: n - layer number
size_t setsim(WSeries< double > &, std::vector< double > *, double=5., double=8., bool saveWF=false)
WSeries< float > v[NIFO_MAX]
void SuperCluster(int ifactor)
double ReadData(double mdcShift, int ifactor)
double threshold(double, double)
param: selected fraction of LTF pixels assuming Gaussian noise param: maximum time delay between dete...