00001
00002
00003
00004
00005
00006
00007
00009
00010 #include <cassert>
00011 #include <csignal>
00012 #include <cmath>
00013 #include <sstream>
00014 #include <algorithm>
00015
00016 #include <Conventions/Munits.h>
00017
00018 #include <EvtDataOverlay/EvtAddRawDigitBlockModule.h>
00019 #include <EvtDataOverlay/EvtOverlaySingleton.h>
00020 #include <Util/WhiteBoardRegistry.h>
00021
00022 #include <MinosObjectMap/MomNavigator.h>
00023 #include <MessageService/MsgService.h>
00024 #include <JobControl/JobCModuleRegistry.h>
00025
00026 #include <Record/SimSnarlHeader.h>
00027 #include <Record/SimSnarlRecord.h>
00028 #include <RawData/RawRecord.h>
00029 #include <RawData/RawDaqSnarlHeader.h>
00030 #include <RawData/RawDigitDataBlock.h>
00031 #include <RawData/RawDigit.h>
00032 #include <RawData/RawVtmTimeInfoBlock.h>
00033
00034 #include <HistMan/HistMan.h>
00035 #include <Calibrator/Calibrator.h>
00036 #include <Plex/PlexHandle.h>
00037
00038 #include <BeamDataUtil/BeamMonSpill.h>
00039
00040 #include "TMath.h"
00041 #include "TString.h"
00042
00043 ClassImp(EvtAddRawDigitBlockModule)
00044
00045
00046
00047 CVSID("$Id: EvtAddRawDigitBlockModule.cxx,v 1.1 2014/06/24 20:25:53 rhatcher Exp $");
00048 JOBMODULE(EvtAddRawDigitBlockModule, "EvtAddRawDigitBlock",
00049 "Adds RawDigitDataBlock to DaqSnarl");
00050
00051
00052
00053
00054 const double tdc_convert = 1000./ 53.1;
00055
00056
00057
00058
00059 EvtAddRawDigitBlockModule::EvtAddRawDigitBlockModule()
00060 : fMaxDeltaT(1.3)
00061 , fSmallT0(20)
00062 {
00063
00064
00065 }
00066
00067
00068
00069 EvtAddRawDigitBlockModule::~EvtAddRawDigitBlockModule()
00070 {
00071 MSG("Hepevt", Msg::kVerbose) << "EvtAddRawDigitBlockModule::Destructor\n";
00072
00073 }
00074
00075
00076
00077 JobCResult EvtAddRawDigitBlockModule::Reco(MomNavigator *mom)
00078 {
00080
00081
00082 RawRecord* rawrec =
00083 dynamic_cast<RawRecord*>(mom->GetFragment("RawRecord",0,"DaqSnarl"));
00084 if ( ! rawrec ) return JobCResult::kFailed;
00085
00086 const RawDaqSnarlHeader* snarlHdr =
00087 dynamic_cast<const RawDaqSnarlHeader*>(rawrec->GetRawHeader());
00088 if ( ! snarlHdr ) return JobCResult::kFailed;
00089
00090 VldTimeStamp vts = snarlHdr->GetVldContext().GetTimeStamp();
00091
00092 EvtOverlaySingleton& eos = EvtOverlaySingleton::Instance();
00093 eos.FindClosestDaqSnarl(vts);
00094 if ( ! eos.IsValidSnarl() ) {
00095 MSG("EvtDataOverlay",Msg::kError) << "invalid DaqSnarl match for "
00096 << vts.AsString("c") << std::endl;
00097 return JobCResult::kFailed;
00098 }
00099 double deltat = eos.GetSnarlTimeDifference();
00100 if ( std::abs(deltat) > fMaxDeltaT ) {
00101 MSG("EvtDataOverlay",Msg::kWarning) << "found closest DaqSnarl record was "
00102 << deltat << " away from "
00103 << vts.AsString("c") << std::endl;
00104 }
00105
00106
00107 const RawDigitDataBlock* origRDDB = eos.GetClosestRawDigitDataBlock();
00108
00109
00110 RawDigitDataBlock* copyRDDB = new RawDigitDataBlock(*origRDDB);
00111
00112
00113
00114
00115 rawrec->AdoptRawBlock(copyRDDB);
00116
00117
00118
00119
00120
00121 eos.FindLowerBoundDaqMonitor(vts);
00122 if ( ! eos.IsValidMonitor() ) {
00123 MSG("EvtDataOverlay",Msg::kError) << "invalid DaqMonitor match for "
00124 << vts.AsString("c") << std::endl;
00125 return JobCResult::kFailed;
00126 }
00127 double deltatm = eos.GetMonitorTimeDifference();
00128
00129 if ( deltatm > fMaxDeltaT || deltatm < 0 ) {
00130 MSG("EvtDataOverlay",Msg::kWarning) << "found lowerbound DaqMonitor record was "
00131 << deltatm << " away from "
00132 << vts.AsString("c") << std::endl;
00133 }
00134
00135
00136 const RawVtmTimeInfoBlock* origRVTIB =
00137 eos.GetLowerBoundRawVtmTimeInfoBlock();
00138 RawVtmTimeInfoBlock* copyRVTIB = new RawVtmTimeInfoBlock(*origRVTIB);
00139 rawrec->AdoptRawBlock(copyRVTIB);
00140
00141 return JobCResult::kAOK;
00142 }
00143
00144
00145
00146 JobCResult EvtAddRawDigitBlockModule::Ana(const MomNavigator *mom)
00147 {
00148
00149
00150 RawRecord* rawrec =
00151 dynamic_cast<RawRecord*>(mom->GetFragment("RawRecord",0,"DaqSnarl"));
00152 if ( ! rawrec ) {
00153 MSG("Overlay",Msg::kWarning)
00154 << "No RawRecord found" << std::endl;
00155 return JobCResult::kFailed;
00156 }
00157
00158 const RawDaqSnarlHeader* snarlHdr =
00159 dynamic_cast<const RawDaqSnarlHeader*>(rawrec->GetRawHeader());
00160 if ( ! snarlHdr ) {
00161 MSG("Overlay",Msg::kWarning)
00162 << "No RawDaqSnarlHeader found" << std::endl;
00163 return JobCResult::kFailed;
00164 }
00165
00166
00167
00168
00169 static bool first = true;
00170 if ( first ) {
00171 int run = snarlHdr->GetRun();
00172 int subrun = snarlHdr->GetSubRun();
00173 BookHist(run,subrun);
00174 first = false;
00175 }
00176
00177 VldContext vldc = snarlHdr->GetVldContext();
00178 Calibrator& cal = Calibrator::Instance();
00179 cal.Reset(vldc);
00180 PlexHandle ph(vldc);
00181
00182 EvtOverlaySingleton& eos = EvtOverlaySingleton::Instance();
00183
00184
00185 HistMan hm("overlay");
00186
00187 std::vector<int> tdc0v, tdc1v;
00188 VldTimeStamp t0_0(VldTimeStamp::GetBOT()), t0_1(VldTimeStamp::GetBOT());
00189
00190 int nrddb = -1;
00191 TIter rbitr = rawrec->GetRawBlockIter();
00192 RawDataBlock* rdb = 0;
00193
00194
00195 while ( (rdb = dynamic_cast<RawDataBlock*>(rbitr())) ) {
00196 RawDigitDataBlock* rddb = dynamic_cast<RawDigitDataBlock*>(rdb);
00197 if ( ! rddb ) continue;
00198
00199 ++nrddb;
00200
00201
00202 const char* name_small = ((nrddb==0) ? "smallT0 block0":"smallT0 block1");
00203 const char* name_big = ((nrddb==0) ? "bigT0 block0" :"bigT0 block1");
00204 TH1D* hsmall = hm.Get<TH1D>(name_small);
00205 TH1D* hbig = hm.Get<TH1D>(name_big);
00206
00207 std::vector<int>& tdcv = ((nrddb==0)? tdc0v : tdc1v);
00208
00209
00210
00211
00212
00213 std::set<VldTimeStamp> t0s = eos.GetCrateT0s(rddb,vldc);
00214 if ( ! t0s.empty() ) {
00215 if (nrddb==0) t0_0 = *(t0s.begin());
00216 else t0_1 = *(t0s.begin());
00217 if ( t0s.size() > 1 ) {
00218 MSG("Overlay",Msg::kWarning)
00219 << "!!!!!!! more than one crateT0 nrddb=" << nrddb << endl;
00220 }
00221 } else {
00222 MSG("Overlay",Msg::kWarning)
00223 << "!!!!!!! no crateT0 nrddb=" << nrddb << endl;
00224 }
00225
00226
00227 TIter rditer = rddb->GetDatumIter();
00228 RawDigit* digit = 0;
00229 int lastCrate = -1;
00230 VldTimeStamp crateT0;
00231 bool isSmallT0 = true;
00232 while ( ( digit = (RawDigit*)rditer() ) ) {
00233 RawChannelId rcid = digit->GetChannel();
00234 ReadoutType::Readout_t digitType = ph.GetReadoutType(rcid);
00235 if ( ReadoutType::kScintStrip != digitType ) continue;
00236
00237
00238 int crate = rcid.GetCrate();
00239 if ( crate != lastCrate ) {
00240 lastCrate = crate;
00241 crateT0 = digit->GetCrateT0();
00242 isSmallT0 = ( abs(crateT0.GetNanoSec()) < fSmallT0 );
00243 }
00244
00245 tdcv.push_back(digit->GetTDC());
00246
00247
00248
00249 double digitTimeNS = cal.GetTimeFromTDC(digit->GetTDC(),rcid)*1.0e9;
00250
00251
00252 digitTimeNS += ( digit->GetCrateT0().GetSec() -
00253 vldc.GetTimeStamp().GetSec() ) * 1.0e9;
00254
00255 double tadj = digitTimeNS -
00256 ( vldc.GetTimeStamp().GetNanoSec() - crateT0.GetNanoSec() );
00257 tadj *= 1.0e-9;
00258
00259 MSG("Overlay",Msg::kDebug)
00260 << "nrddb " << nrddb << "|" << ((isSmallT0)?"s":"B")
00261 << " tadj = " << setw(8) << tadj/Munits::microsecond << " usec "
00262 << " TDC " << digit->GetTDC()
00263
00264
00265
00266 << " " << rcid.AsString("c")
00267 << endl;
00268
00269 if ( isSmallT0 ) hsmall->Fill(tadj);
00270 else hbig->Fill(tadj);
00271
00272 }
00273
00274 if ( lastCrate != -1 ) {
00275
00276 TH1I* hseen = hm.Get<TH1I>("rddb seen");
00277 if ( isSmallT0 ) hseen->Fill(name_small,1.0);
00278 else hseen->Fill(name_big,1.0);
00279
00280
00281
00282
00283
00284
00285
00286
00287
00288
00289
00290
00291
00292
00293
00294 }
00295
00296 }
00297
00298
00299
00300 if ( t0_0 != t0_1 ) {
00301 MSG("Overlay",Msg::kWarning)
00302 << "!!!!!!! RDDB blocks had different crateT0 "
00303 << t0_0.AsString("c") << " " << t0_1.AsString("c") << endl;
00304 }
00305
00306 std::sort(tdc0v.begin(),tdc0v.end());
00307 std::sort(tdc1v.begin(),tdc1v.end());
00308
00309 int maxtdc = *(tdc1v.rbegin());
00310 int mintdc = *(tdc1v.begin());
00311
00312 int snarl = snarlHdr->GetSnarl();
00313 std::ostringstream s;
00314 s << "snarl" << snarl << "_";
00315 std::string name0 = s.str() + "0";
00316 std::string name1 = s.str() + "1";
00317
00318 MSG("Overlay",Msg::kSynopsis)
00319 << "Snarl " << snarl << " TDC min/max " << mintdc << " " << maxtdc
00320 << " range " << (maxtdc-mintdc) << endl;
00321
00322 #undef HIST_INDIVIDUAL
00323 #ifdef HIST_INDIVIDUAL
00324 HistMan hmMany("individual");
00325 TH1I* htdc0 = hmMany.Book<TH1I>(name0.c_str(),name0.c_str(),(maxtdc-mintdc)+1,mintdc,maxtdc+1);
00326 TH1I* htdc1 = hmMany.Book<TH1I>(name1.c_str(),name1.c_str(),(maxtdc-mintdc)+1,mintdc,maxtdc+1);
00327
00328 for (size_t i=0; i<tdc0v.size(); ++i) htdc0->Fill(tdc0v[i]);
00329 for (size_t i=0; i<tdc1v.size(); ++i) htdc1->Fill(tdc1v[i]);
00330 #endif
00331
00332
00333 int index = 0;
00334 minos::GetFromWhiteBoard("EvtInfo","TimeIndex",index);
00335 double pots = eos.GetClosestSpillPOTs(true);
00336 const BeamMonSpill* bms = eos.GetClosestBeamMonSpill();
00337
00338 HistMan hmPOTs(fHM_pots_dir.c_str());
00339 hmPOTs.Get<TH1D>("POTs")->Fill(index,pots);
00340 TH2D* hbi = hmPOTs.Get<TH2D>("BatchIntensity");
00341 if ( tdc0v.size() > 0 ) {
00342 double x = 100. + log((double)(tdc0v.size()));
00343 hbi->Fill(index,7,x);
00344 }
00345 if ( bms ) {
00346 for (int ib = 0; ib < 6; ++ib) {
00347 double batchi = bms->fBpmInt[ib];
00348 hbi->Fill(index,ib,batchi);
00349 }
00350 }
00351
00352 return JobCResult::kPassed;
00353 }
00354
00355
00356 const Registry& EvtAddRawDigitBlockModule::DefaultConfig() const
00357 {
00358
00359
00360
00361 static Registry r;
00362
00363
00364 std::string name = this->GetName();
00365 name += ".config.default";
00366 r.SetName(name.c_str());
00367
00368
00369 r.UnLockValues();
00370
00371 r.Set("MaxDeltaT",1.3);
00372
00373 r.LockValues();
00374
00375 return r;
00376 }
00377
00378
00379
00380 void EvtAddRawDigitBlockModule::Config(const Registry& r)
00381 {
00382
00383
00384
00385
00386 double tmpd;
00387
00388
00389 if (r.Get("MaxDeltaT", tmpd)) { fMaxDeltaT = tmpd; }
00390
00391 }
00392
00393
00394
00395 void EvtAddRawDigitBlockModule::BeginJob()
00396 {
00397
00398 }
00399
00400 void EvtAddRawDigitBlockModule::BookHist(int run, int subrun)
00401 {
00402
00403 static bool first = true;
00404 if ( ! first ) return;
00405 first = false;
00406
00407
00408 HistMan hm("overlay");
00409
00410 double upper = 20000. * 1.0e-9;
00411 int ngang = 1;
00412 int nbins = TMath::Nint((double)upper/ (tdc_convert*1.0e-9)) / ngang;
00413
00414 TH1I* hseen = hm.Book<TH1I>("rddb seen","rddb seen",4,0,4);
00415 hseen->SetBit(TH1::kCanRebin);
00416
00417 const char* hnames[4] = { "smallT0 block0", "smallT0 block1",
00418 "bigT0 block0", "bigT0 block1" };
00419
00420 for (int ihn=0; ihn<4; ++ihn) {
00421 const char* hname = hnames[ihn];
00422 std::string htitle(hname);
00423 htitle += "; (sec); # digits";
00424 hm.Book<TH1D>(hname,htitle.c_str(),nbins,0,upper);
00425
00426 hseen->Fill(hname,0.0);
00427 }
00428 hseen->LabelsOption("a");
00429
00430 MSG("Overlay",Msg::kInfo)
00431 <<"EvtAddRawDigitBlockModule BOOK block[0|1] [small|big]T0 "
00432 << "(" << nbins << ",0," << upper <<")" << endl;
00433
00434
00435
00436 fHM_pots_dir = Form("r%ds%d POTs",run,subrun);
00437 HistMan hmPOTs(fHM_pots_dir.c_str());
00438
00439
00440
00441 EvtOverlaySingleton& eos = EvtOverlaySingleton::Instance();
00442 size_t ntimes = eos.NumUserTimes();
00443
00444 hmPOTs.Book<TH1D>("POTs","POTs",ntimes,0,ntimes);
00445 hmPOTs.Book<TH2D>("BatchIntensity","BatchIntensity",ntimes,0,ntimes,8,0,8);
00446
00447 }
00448
00449
00450
00451 void EvtAddRawDigitBlockModule::EndJob()
00452 {
00453
00454 MSG("Hepevt",Msg::kInfo)
00455 << "EvtAddRawDigitBlockModule Added " << endl;
00456
00457 }
00458
00459