00001
00002
00003
00004
00005
00006
00007
00008 #include "TH1.h"
00009 #include "TH1D.h"
00010 #include "TH2D.h"
00011 #include "TProfile2D.h"
00012 #include "TMath.h"
00013 #include "TDirectory.h"
00014 #include "TROOT.h"
00015 #include <iomanip>
00016 #include "TChain.h"
00017 #include "TChainElement.h"
00018 #include "TLeaf.h"
00019 #include <math.h>
00020 #include <fstream>
00021 #include <sstream>
00022 #include <string>
00023 #include <iostream>
00024 #include <iomanip>
00025
00026 #include "Conventions/Detector.h"
00027 #include "Conventions/SimFlag.h"
00028 #include "MessageService/MsgService.h"
00029 #include "NtupleUtils/NuCutter.h"
00030 #include "NtupleUtils/NuEvent.h"
00031 #include "NtupleUtils/NuInputEvents.h"
00032 #include "NtupleUtils/NuLibrary.h"
00033 #include "NtupleUtils/NuSystematic.h"
00034 #include "NuMuBar/NuTransition.h"
00035 #include "NtupleUtils/NuXMLConfig.h"
00036 #include "NuMuBar/DataStabilityAnalysisClass.h"
00037 #include "Validity/VldTimeStamp.h"
00038 #include "BeamDataUtil/BDSpillAccessor.h"
00039 #include "BeamDataUtil/BeamMonSpill.h"
00040 #include "Util/MsgUtil.h"
00041
00042 using std::vector;
00043
00044 CVSID("$Id: DataStabilityAnalysisClass.cxx,v 1.6 2015/09/30 15:02:09 ashley90 Exp $");
00045
00047 DataStabilityAnalysisClass::DataStabilityAnalysisClass(NuXMLConfig* xmlConfig)
00048 {
00049 if (!xmlConfig){
00050 MSG("DataStabilityAnalysisClass",Msg::kWarning)
00051 << "Analysis object configured with a non-existent xml object."
00052 << endl
00053 << "Asserting now."
00054 << endl;
00055 assert(false);
00056 }
00057 fxmlConfig = xmlConfig;
00058
00059
00060 TString anaVer = fxmlConfig->AnaVersionString();
00061
00062 fCutter = new NuCutter(anaVer);
00063 timeTree = new TTree("timeTree","Tree of times");
00064
00065
00066 fAnaVersion = static_cast<NuCuts::NuAnaVersion_t>(fCutter->AnaVersion());
00067
00068
00069 fsyst = new NuSystematic(*xmlConfig);
00070 fsyst->SetNuBarSelector(fCutter->GetCut());
00071 fsyst->SetCCSelector(fCutter->GetCut());
00072 fsyst->SetNCSelector(fCutter->GetCut());
00073
00074 fLib = &NuLibrary::Instance();
00075
00076
00077 prevSnarl = 0.0;
00078 tmpSubRun = 0;
00079 tmprun = 0;
00080 newSubRun = false;
00081
00082
00083
00084
00085 FixedTime = 0;
00086
00087
00088
00089
00090
00091
00092 StartTime = 1062547200;
00093 EndTime = 1378080000;
00094
00095 fHistogramsCreated = false;
00096 fTreeCreated = false;
00097
00098 DSTpotCounter = 0.0;
00099 DSTspillcounter = 0;
00100 delta = 0;
00101
00102
00103 tselectionevent = -999;
00104 ttimeRelToSpill = -999;
00105 ttimeSec = -999;
00106 ttimeNanoSec = -999;
00107 tnearestSpillSec = -999;
00108 tnearestSpillNanosec = -999;
00109 tnshw = -999;
00110 trun = -999;
00111 tsubRun = -999;
00112 tsnarl = -999;
00113 ttrkEn = -999.;
00114 ttrkEnRange = -999.;
00115 ttrkEnCurv = -999.;
00116 tshwEn = -999.;
00117 tshwEnLinCCCor = -999.;
00118 tshwEnkNN = -999.;
00119 txEvtVtx = -999.;
00120 tyEvtVtx = -999.;
00121 tzEvtVtx = -999.;
00122 fTor101 = -999.;
00123 fTr101d = -999.;
00124 fTortgt = -999.;
00125 fTrtgtd = -999.;
00126 fHornCur = -999.;
00127 troID = -999.;
00128 fTargProfX = -999.;
00129 fTargProfY = -999.;
00130 fProfWidX = -999.;
00131 fProfWidY = -999.;
00132 fHadInt = 0.;
00133 fMuInt1 = 0.;
00134 fMuInt2 = 0.;
00135 fMuInt3 = 0.;
00136 totBpm = 0.;
00137 whichBatch = -1;
00138 numBatches = 0;
00139
00140 beamTimeSec = 0.;
00141 beamTimeNSec = 0.;
00142
00143 batchBin1summer = 0.0;
00144 }
00145
00147 DataStabilityAnalysisClass::~DataStabilityAnalysisClass()
00148 {
00149 if (fCutter) {
00150 delete fCutter;
00151 fCutter = 0;
00152 }
00153 fLib = 0;
00154 this->DeleteHistograms();
00155 }
00157 NuInputEvents& DataStabilityAnalysisClass::GetEventListForLooping()
00158 {
00159 NuInputEvents& input = this->DoIO(0, "", "null");
00160 fxmlConfig->Write();
00161 input.ResetNuEventLoopPositionToStart();
00162 return input;
00163 }
00164
00166 void DataStabilityAnalysisClass::MakePlots()
00167 {
00168 NuInputEvents& input = this->GetEventListForLooping();
00169
00170 BDSpillAccessor& bsAccess = BDSpillAccessor::Get();
00171 BeamMonSpill *beamSpill = 0x0;
00172
00173 for (Int_t i = 0; i < input.GetEntriesNuEvent(); ++i){
00174
00175 NuUtilities::ProgressBar(i, input.GetEntriesNuEvent(), 5);
00176 NuEvent &nu = const_cast<NuEvent&>(input.GetNextNuEvent(Msg::kDebug));
00177
00178 if (!fHistogramsCreated){ this->CreateHistograms(nu); }
00179 if (!fTreeCreated){ this->CreateTree(nu); }
00180
00181 this->PrepareEventForPlotting(nu);
00182
00183
00184
00185
00186 if(nu.run != tmprun){
00187 FixedTime = nu.timeSeconds;
00188 }
00189
00190
00191
00192
00193
00194
00195
00196 if(nu.run != tmprun){
00197 tmprun = nu.run;
00198 }
00199
00200
00201 fCutter->MakeCuts(nu);
00202
00203
00204
00205 Bool_t passedSelectionCut = false;
00206
00207
00208 if( fCutter->AnaVersion() == NuCuts::kNCPRLNC ||
00209 fCutter->AnaVersion() == NuCuts::kCC1070Std ||
00210 fCutter->AnaVersion() == NuCuts::kCC0720Std ){
00211
00212 if(fCutter->Passed() == true){
00213 passedSelectionCut = true;
00214 }
00215 else{
00216 passedSelectionCut = false;
00217 }
00218 }
00219 if( fCutter->AnaVersion() == NuCuts:: kNMB0720Bravo ){
00220 if(fCutter->Passed() == true && nu.charge == +1){
00221 passedSelectionCut = true;
00222 }
00223 else{
00224 passedSelectionCut = false;
00225 }
00226 }
00227
00228
00229
00230
00231
00232 if ( fCutter->Passed() ){
00233
00234
00235 if( fCutter->AnaVersion() == NuCuts::kNCPRLNC ||
00236 fCutter->AnaVersion() == NuCuts::kCC1070Std ||
00237 fCutter->AnaVersion() == NuCuts::kCC0720Std ){
00238
00239
00240
00241
00242
00243 this->FillPlots(nu);
00244 }
00245
00246
00247 if( fCutter->AnaVersion() == NuCuts:: kNMB0720Bravo ){
00248
00249
00250 if(nu.charge == +1){
00251 this->FillPlots(nu);
00252
00253 }
00254 }
00255 }
00256 if (fCutter->Failed()){
00257 }
00258 }
00259
00260 this->WritePlots();
00261 fCutter->PrintSummary();
00262 fCutter->PrintNMinusOneSummary();
00263 fOutFile->Close();
00264 }
00265
00267 void DataStabilityAnalysisClass::PrepareEventForPlotting(NuEvent& event) const
00268 {
00269 event.anaVersion = fAnaVersion;
00270 fxmlConfig->ConfigureWeights(event);
00271
00272
00273
00274
00275
00276
00277
00278 if (event.simFlag!=SimFlag::kData) {
00279 fsyst->Shift(event);
00280
00281 }
00282
00283
00284 if(NuCuts::kNCPRLNC == fAnaVersion){
00285
00286 }
00287
00288 }
00290 void DataStabilityAnalysisClass::CreateTree(const NuEvent& event)
00291 {
00292 std::cout << "Creating TTree." << std::endl;
00293
00294
00295 Util::MsgLevel("Dbi","Fatal");
00296 Util::MsgLevel("BDU","Fatal");
00297
00298 if (event.simFlag == SimFlag::kMC){
00299
00300 }
00301
00302
00303 timeTree->Branch("nshw", &tnshw, "nshw/I");
00304 timeTree->Branch("selectionevent", &tselectionevent, "selectionevent/I");
00305 timeTree->Branch("timeRelToSpill",&ttimeRelToSpill,"timeRelToSpill/D");
00306 timeTree->Branch("timeSec",&ttimeSec,"timeSec/I");
00307 timeTree->Branch("timeNanoSec",&ttimeNanoSec,"timeNanoSec/I");
00308 timeTree->Branch("nearestSpillSec",&tnearestSpillSec,"nearestSpillSec/I");
00309 timeTree->Branch("nearestSpillNanosec",&tnearestSpillNanosec,"nearestSpillNanosec/I");;
00310 timeTree->Branch("run",&trun,"run/I");
00311 timeTree->Branch("subRun",&tsubRun,"subRun/I");
00312 timeTree->Branch("snarl",&tsnarl,"snarl/I");
00313 timeTree->Branch("trkEn",&ttrkEn,"trkEn/F");
00314 timeTree->Branch("trkEnRange",&ttrkEnRange,"trkEnRange/F");
00315 timeTree->Branch("trkEnCurv",&ttrkEnCurv,"trkEnCurv/F");
00316 timeTree->Branch("shwEn",&tshwEn,"shwEn/F");
00317 timeTree->Branch("shwEnLinCCCor",&tshwEnLinCCCor,"shwEnLinCCCor/F");
00318 timeTree->Branch("shwEnkNN",&tshwEnkNN,"shwEnkNN/F");
00319 timeTree->Branch("roID", &troID, "roID/F");
00320 timeTree->Branch("xEvtVtx",&txEvtVtx,"xEvtVtx/F");
00321 timeTree->Branch("yEvtVtx",&tyEvtVtx,"yEvtVtx/F");
00322 timeTree->Branch("zEvtVtx",&tzEvtVtx,"zEvtVtx/F");
00323
00324
00325 timeTree->Branch("batchPot", batchPot,"batchPot[6]/F");
00326 timeTree->Branch("numBatches",&numBatches,"numBatches/I");
00327 timeTree->Branch("whichBatch",&whichBatch,"whichBatch/I");
00328 timeTree->Branch("spillTimeSec",&beamTimeSec,"spillTimeSec/l");
00329 timeTree->Branch("spillTimeNSec",&beamTimeNSec,"spillTimeNSec/l");
00330 timeTree->Branch("fTor101",&fTor101,"fTor101/F");
00331 timeTree->Branch("fTr101d",&fTr101d,"fTr101d/F");
00332 timeTree->Branch("fTrtgtd",&fTrtgtd,"fTrtgtd/F");
00333 timeTree->Branch("fTortgt",&fTortgt,"fTortgt/F");
00334 timeTree->Branch("fHornCur",&fHornCur,"fHornCur/F");
00335 timeTree->Branch("fTargBpmX",fTargBpmX,"fTargBpmX[6]/F");
00336 timeTree->Branch("fTargBpmY",fTargBpmY,"fTargBpmY[6]/F");
00337 timeTree->Branch("fBpmInt",fBpmInt,"fBpmInt[6]/F");
00338 timeTree->Branch("fTargProfX",&fTargProfX,"fTargProfX/F");
00339 timeTree->Branch("fTargProfY",&fTargProfY,"fTargProfY/F");
00340 timeTree->Branch("fProfWidX",&fProfWidX,"fProfWidX/F");
00341 timeTree->Branch("fProfWidY",&fProfWidY,"fProfWidY/F");
00342 timeTree->Branch("fHadInt",&fHadInt,"fHadInt/F");
00343 timeTree->Branch("fMuInt1",&fMuInt1,"fMuInt1/F");
00344 timeTree->Branch("fMuInt2",&fMuInt2,"fMuInt2/F");
00345 timeTree->Branch("fMuInt3",&fMuInt3,"fMuInt3/F");
00346
00347 fTreeCreated = true;
00348 }
00350 void DataStabilityAnalysisClass::CreateHistograms(const NuEvent& event)
00351 {
00352 std::cout << "Creating histograms" << std::endl;
00353
00354
00355 const vector<Double_t> recoBins = NuUtilities::RecoBins(NuBinningScheme::kSterile);
00356 const Int_t numRecoBins = recoBins.size() - 1;
00357 Double_t* recoBinsArray = new Double_t[numRecoBins + 1];
00358 Int_t i = 0;
00359
00360
00361 for (vector<Double_t>::const_iterator itBin = recoBins.begin();
00362 itBin != recoBins.end();
00363 ++itBin, ++i){
00364
00365 recoBinsArray[i] = *itBin;
00366 }
00367
00368 hRecoEnergyVsTime = new TH2D("hRecoEnergyVsTime","",3652, StartTime, EndTime, 400, 0, 200);
00369
00370
00371
00372 hRecoEnergyAll = new TH1D("hRecoEnergyAll","",numRecoBins, recoBinsArray);
00373 hRecoEnergyNQ = new TH1D("hRecoEnergyNQ","", 400, 0, 200);
00374 hRecoEnergyPQ = new TH1D("hRecoEnergyPQ","", 400, 0, 200);
00375
00376 hEventsPerTimePerPOTLT3GeV = new TH1D("hEventsPerTimePerPOTLT3GeV", "", 3652, StartTime, EndTime);
00377 hEventsPerTimePerPOTGT3GeV = new TH1D("hEventsPerTimePerPOTGT3GeV", "", 3652, StartTime, EndTime);
00378 hEventsPerTimePerPOTLT6GeV = new TH1D("hEventsPerTimePerPOTLT6GeV", "", 3652, StartTime, EndTime);
00379 hEventsPerTimePerPOTGT6GeV = new TH1D("hEventsPerTimePerPOTGT6GeV", "", 3652, StartTime, EndTime);
00380 hEventsPerTimePerPOTLT8GeV = new TH1D("hEventsPerTimePerPOTLT8GeV", "", 3652, StartTime, EndTime);
00381 hEventsPerTimePerPOTGT8GeV = new TH1D("hEventsPerTimePerPOTGT8GeV", "", 3652, StartTime, EndTime);
00382 hEventsPerTimePerPOTLT14GeV = new TH1D("hEventsPerTimePerPOTLT14GeV", "", 3652, StartTime, EndTime);
00383 hEventsPerTimePerPOTGT14GeV = new TH1D("hEventsPerTimePerPOTGT14GeV", "", 3652, StartTime, EndTime);
00384
00385
00386 hRecoEnergyAllmuon = new TH1D("hRecoEnergyAllmuon","",400,0,200);
00387 hEventsPerTimePerPOTLT3GeVmuon = new TH1D("hEventsPerTimePerPOTLT3GeVmuon", "", 3652, StartTime, EndTime);
00388 hEventsPerTimePerPOTGT3GeVmuon = new TH1D("hEventsPerTimePerPOTGT3GeVmuon", "", 3652, StartTime, EndTime);
00389 hEventsPerTimePerPOTLT6GeVmuon = new TH1D("hEventsPerTimePerPOTLT6GeVmuon", "", 3652, StartTime, EndTime);
00390 hEventsPerTimePerPOTGT6GeVmuon = new TH1D("hEventsPerTimePerPOTGT6GeVmuon", "", 3652, StartTime, EndTime);
00391 hEventsPerTimePerPOTLT8GeVmuon = new TH1D("hEventsPerTimePerPOTLT8GeVmuon", "", 3652, StartTime, EndTime);
00392 hEventsPerTimePerPOTGT8GeVmuon = new TH1D("hEventsPerTimePerPOTGT8GeVmuon", "", 3652, StartTime, EndTime);
00393 hEventsPerTimePerPOTLT14GeVmuon = new TH1D("hEventsPerTimePerPOTLT14GeVmuon", "", 3652, StartTime, EndTime);
00394 hEventsPerTimePerPOTGT14GeVmuon = new TH1D("hEventsPerTimePerPOTGT14GeVmuon", "", 3652, StartTime, EndTime);
00395
00396
00397 hbatchPOT_1batchmode = new TH2D( "batchPOT_1batchmode", "batchPOT_1batchmode; Batch Index; Batch Intensity (x10^{12} POT)", 7, 0, 7, 200, 0, 20);
00398 hbatchPOT_2batchmode = new TH2D( "batchPOT_2batchmode", "batchPOT_2batchmode; Batch Index; Batch Intensity (x10^{12} POT)", 7, 0, 7, 200, 0, 20);
00399 hbatchPOT_3batchmode = new TH2D( "batchPOT_3batchmode", "batchPOT_3batchmode; Batch Index; Batch Intensity (x10^{12} POT)", 7, 0, 7, 200, 0, 20);
00400 hbatchPOT_4batchmode = new TH2D( "batchPOT_4batchmode", "batchPOT_4batchmode; Batch Index; Batch Intensity (x10^{12} POT)", 7, 0, 7, 200, 0, 20);
00401 hbatchPOT_5batchmode = new TH2D( "batchPOT_5batchmode", "batchPOT_5batchmode; Batch Index; Batch Intensity (x10^{12} POT)", 7, 0, 7, 200, 0, 20);
00402 hbatchPOT_6batchmode = new TH2D( "batchPOT_6batchmode", "batchPOT_6batchmode; Batch Index; Batch Intensity (x10^{12} POT)", 7, 0, 7, 200, 0, 20);
00403
00404 hBatchPOTSum_1batchmode = new TH1D( "hBatchPOTSum_1batchmode","hBatchPOTSum_1batchmode; Batch Index; Batch POT (x10^{12} POT) ", 7, 0, 7);
00405 hBatchPOTSum_2batchmode = new TH1D( "hBatchPOTSum_2batchmode","hBatchPOTSum_2batchmode; Batch Index; Batch POT (x10^{12} POT) ", 7, 0, 7);
00406 hBatchPOTSum_3batchmode = new TH1D( "hBatchPOTSum_3batchmode","hBatchPOTSum_3batchmode; Batch Index; Batch POT (x10^{12} POT) ", 7, 0, 7);
00407 hBatchPOTSum_4batchmode = new TH1D( "hBatchPOTSum_4batchmode","hBatchPOTSum_4batchmode; Batch Index; Batch POT (x10^{12} POT) ", 7, 0, 7);
00408 hBatchPOTSum_5batchmode = new TH1D( "hBatchPOTSum_5batchmode","hBatchPOTSum_5batchmode; Batch Index; Batch POT (x10^{12} POT) ", 7, 0, 7);
00409 hBatchPOTSum_6batchmode = new TH1D( "hBatchPOTSum_6batchmode","hBatchPOTSum_6batchmode; Batch Index; Batch POT (x10^{12} POT) ", 7, 0, 7);
00410
00411
00412 hSelectionBatchesAll = new TH1D("hSelectionBatchesAll" , "", 4000, 0, 12e-6);
00413 hSelectionBatches_BatchMode6 = new TH1D("hSelectionBatches_BatchMode6", "", 4000, 0, 12e-6);
00414 hSelectionBatches_BatchMode5 = new TH1D("hSelectionBatches_BatchMode5", "", 4000, 0, 12e-6);
00415 hSelectionBatches_BatchMode4 = new TH1D("hSelectionBatches_BatchMode4", "", 4000, 0, 12e-6);
00416 hSelectionBatches_BatchMode3 = new TH1D("hSelectionBatches_BatchMode3", "", 4000, 0, 12e-6);
00417 hSelectionBatches_BatchMode2 = new TH1D("hSelectionBatches_BatchMode2", "", 4000, 0, 12e-6);
00418 hSelectionBatches_BatchMode1 = new TH1D("hSelectionBatches_BatchMode1", "", 4000, 0, 12e-6);
00419
00420
00421 hCalibrateBatchPOT = new TH2D("hCalibrateBatchPOT","CalibrateBatchPOT; fTrtgtd (x10^{12} POT); #Sigma fBpmInt", 250, 0, 40, 500, 0, 4000);
00422
00423
00424 hVertexZvsEnergy = new TH2D("hVertexZvsEnergy","",400,0,200,100,0,8);
00425
00426
00427 hPotPerIntensity = new TH1D("hPotPerIntensity", "hPotPerIntensity; Intensity per Spill (x10^{12} POT); Summed POT x10^{12}", 450, 0, 50);
00428
00429
00430 hPotPerIntensity_5BatchMode_Batch[0] = new TH1D("hPotPerIntensity_5BatchMode_Batch1",
00431 "hPotPerIntensity_5BatchMode_Batch1; Intensity per Batch (x10^{12} POT); Summed POT x10^{12}", 400, 0, 20);
00432 hPotPerIntensity_5BatchMode_Batch[1] = new TH1D("hPotPerIntensity_5BatchMode_Batch2",
00433 "hPotPerIntensity_5BatchMode_Batch2; Intensity per Batch (x10^{12} POT); Summed POT x10^{12}", 400, 0, 20);
00434 hPotPerIntensity_5BatchMode_Batch[2] = new TH1D("hPotPerIntensity_5BatchMode_Batch3",
00435 "hPotPerIntensity_5BatchMode_Batch3; Intensity per Batch (x10^{12} POT); Summed POT x10^{12}", 400, 0, 20);
00436 hPotPerIntensity_5BatchMode_Batch[3] = new TH1D("hPotPerIntensity_5BatchMode_Batch4",
00437 "hPotPerIntensity_5BatchMode_Batch4; Intensity per Batch (x10^{12} POT); Summed POT x10^{12}", 400, 0, 20);
00438 hPotPerIntensity_5BatchMode_Batch[4] = new TH1D("hPotPerIntensity_5BatchMode_Batch5",
00439 "hPotPerIntensity_5BatchMode_Batch5; Intensity per Batch (x10^{12} POT); Summed POT x10^{12}", 400, 0, 20);
00440
00441 hPotPerIntensity_6BatchMode_Batch[0] = new TH1D("hPotPerIntensity_6BatchMode_Batch1",
00442 "hPotPerIntensity_6BatchMode_Batch1; Intensity per Batch (x10^{12} POT); Summed POT x10^{12}", 400, 0, 20);
00443 hPotPerIntensity_6BatchMode_Batch[1] = new TH1D("hPotPerIntensity_6BatchMode_Batch2",
00444 "hPotPerIntensity_6BatchMode_Batch2; Intensity per Batch (x10^{12} POT); Summed POT x10^{12}", 400, 0, 20);
00445 hPotPerIntensity_6BatchMode_Batch[2] = new TH1D("hPotPerIntensity_6BatchMode_Batch3",
00446 "hPotPerIntensity_6BatchMode_Batch3; Intensity per Batch (x10^{12} POT); Summed POT x10^{12}", 400, 0, 20);
00447 hPotPerIntensity_6BatchMode_Batch[3] = new TH1D("hPotPerIntensity_6BatchMode_Batch4",
00448 "hPotPerIntensity_6BatchMode_Batch4; Intensity per Batch (x10^{12} POT); Summed POT x10^{12}", 400, 0, 20);
00449 hPotPerIntensity_6BatchMode_Batch[4] = new TH1D("hPotPerIntensity_6BatchMode_Batch5",
00450 "hPotPerIntensity_6BatchMode_Batch5; Intensity per Batch (x10^{12} POT); Summed POT x10^{12}", 400, 0, 20);
00451 hPotPerIntensity_6BatchMode_Batch[5] = new TH1D("hPotPerIntensity_6BatchMode_Batch6",
00452 "hPotPerIntensity_6BatchMode_Batch6; Intensity per Batch (x10^{12} POT); Summed POT x10^{12}", 400, 0, 20);
00453
00454
00455
00456 IntensityPerTime = new TH2D("IntensityPerTime", "IntensityPerTime; Time; Spill Intensity (x10^{12} ", 3652, StartTime, EndTime, 450, 0, 50);
00457 hNumberofGoodSpillsPerIntensity = new TH1D("NumberofGoodSpillsPerIntensity", "NumberofGoodSpillsPerIntensity; Intensity per Spill (x10^{12} POT); # Spills", 400, 0, 45);
00458
00459 IntensityPerTime_5BatchMode_Batch[0] = new TH2D("IntensityPerTime_5BatchMode_Batch1",
00460 "IntensityPerTime_5BatchMode_Batch1; Time; Spill Intensity (x10^{12} ", 3652, StartTime, EndTime, 400, 0, 20);
00461 IntensityPerTime_5BatchMode_Batch[1] = new TH2D("IntensityPerTime_5BatchMode_Batch2",
00462 "IntensityPerTime_5BatchMode_Batch2; Time; Spill Intensity (x10^{12} ", 3652, StartTime, EndTime, 400, 0, 20);
00463 IntensityPerTime_5BatchMode_Batch[2] = new TH2D("IntensityPerTime_5BatchMode_Batch3",
00464 "IntensityPerTime_5BatchMode_Batch3; Time; Spill Intensity (x10^{12} ", 3652, StartTime, EndTime, 400, 0, 20);
00465 IntensityPerTime_5BatchMode_Batch[3] = new TH2D("IntensityPerTime_5BatchMode_Batch4",
00466 "IntensityPerTime_5BatchMode_Batch4; Time; Spill Intensity (x10^{12} ", 3652, StartTime, EndTime, 400, 0, 20);
00467 IntensityPerTime_5BatchMode_Batch[4] = new TH2D("IntensityPerTime_5BatchMode_Batch5",
00468 "IntensityPerTime_5BatchMode_Batch5; Time; Spill Intensity (x10^{12} ", 3652, StartTime, EndTime, 400, 0, 20);
00469
00470 IntensityPerTime_6BatchMode_Batch[0] = new TH2D("IntensityPerTime_6BatchMode_Batch1",
00471 "IntensityPerTime_6BatchMode_Batch1; Time; Spill Intensity (x10^{12} ", 3652, StartTime, EndTime, 400, 0, 20);
00472 IntensityPerTime_6BatchMode_Batch[1] = new TH2D("IntensityPerTime_6BatchMode_Batch2",
00473 "IntensityPerTime_6BatchMode_Batch2; Time; Spill Intensity (x10^{12} ", 3652, StartTime, EndTime, 400, 0, 20);
00474 IntensityPerTime_6BatchMode_Batch[2] = new TH2D("IntensityPerTime_6BatchMode_Batch3",
00475 "IntensityPerTime_6BatchMode_Batch3; Time; Spill Intensity (x10^{12} ", 3652, StartTime, EndTime, 400, 0, 20);
00476 IntensityPerTime_6BatchMode_Batch[3] = new TH2D("IntensityPerTime_6BatchMode_Batch4",
00477 "IntensityPerTime_6BatchMode_Batch4; Time; Spill Intensity (x10^{12} ", 3652, StartTime, EndTime, 400, 0, 20);
00478 IntensityPerTime_6BatchMode_Batch[4] = new TH2D("IntensityPerTime_6BatchMode_Batch5",
00479 "IntensityPerTime_6BatchMode_Batch5; Time; Spill Intensity (x10^{12} ", 3652, StartTime, EndTime, 400, 0, 20);
00480 IntensityPerTime_6BatchMode_Batch[5] = new TH2D("IntensityPerTime_6BatchMode_Batch6",
00481 "hPotPerIntensity_6BatchMode_Batch6; Time; Spill Intensity (x10^{12} ", 3652, StartTime, EndTime, 400, 0, 20);
00482
00483
00484 hedgeActivityPH = new TH1D("hedgeActivityPH","",30,0,10000);
00485 hedgeActivityStrips = new TH1D("hedgeActivityStrips","",20,0,20);
00486 hminTimeSeparation = new TH1D("hminTimeSeparation","",100,-1000,1000);
00487 hoppEdgePH = new TH1D("hoppEdgePH","",30,0,5000);
00488 hoppEdgeStrips = new TH1D("hoppEdgeStrips","",20,0,16);
00489 hcloseTimeDeltaZ = new TH1D("hcloseTimeDeltaZ","",100,-10,10);
00490 hnstripEvt = new TH1D("hnstripEvt","",100,0,300);
00491 hstripsOverSqPlanes = new TH1D("hstripsOverSqPlanes","",80,0,0.4);
00492
00493
00494 hydistribution = new TH1D("hydistribution","",200,0,1);
00495 hplanesEvent = new TH1D("hplanesEvent","",200,0,200);
00496 hPhEvent = new TH1D("hPhEvent","",300,0,6e5);
00497 hPhPerPlane = new TH1D("hPhPerPlane","",200,0,1e4);
00498 hPhPerStrip = new TH1D("hPhPerStrip","",200,0,1e4);
00499 hStripsPerPlane = new TH1D("hStripsPerPlane","",15,0,15);
00500 hTracksInEvent = new TH1D("hTracksInEvent","",8,0,8);
00501 hShowersInEvent = new TH1D("hShowersInEvent","",10,0,10);
00502 hEventXVertex = new TH1D("hEventXVertex","",100,-1,3);
00503 hEventYVertex = new TH1D("hEventYVertex","",100,-2,2);
00504 hEventZVertex = new TH1D("hEventZVertex","",100,0,8);
00505 hEventRVertex = new TH1D("hEventRVertex","",100,0,3);
00506
00507
00508 hTrkLength = new TH1D("hTrkLength","",280,0,280);
00509 hPlanesCrossedInTrk = new TH1D("hPlanesCrossedInTrk","",300,0,300);
00510 hTrkCosTheta = new TH1D("hTrkCosTheta","",100,-1,1);
00511 hTrkdzds = new TH1D("hTrkdzds","",100,0,1);
00512 hTrkPh = new TH1D("hTrkPh","",50,0,20e4);
00513 hTrkMomentumRange = new TH1D("hTrkMomentumRange","",400,0,100);
00514 hTrkPhPerPlane = new TH1D("hTrkPhPerPlane","",20,0,2);
00515 hTrkqpSigma = new TH1D("hTrkqpSigma","",400,0,100);
00516 hTrkPassedFit = new TH1D("hTrkPassedFit","",2,0,2);
00517 hTrkChi2PerNdf = new TH1D("hTrkChi2PerNdf","",100,0,6);
00518 hTrkXVertex = new TH1D("hTrkXVertex","",50,-1,4);
00519 hTrkYVertex = new TH1D("hTrkYVertex","",50,-2,3);
00520 hTrkZVertex = new TH1D("hTrkZVertex","",200,0,8);
00521
00522
00523 hPlanesCrossedInShower = new TH1D("hPlanesCrossedInShower","",100,0,100);
00524 hShowerEnergy = new TH1D("hShowerEnergy","",100,0,50);
00525 hXShwVtx = new TH1D("hXShwVtx","",50,-1,4);
00526 hYShwVtx = new TH1D("hYShwVtx","",50,-2,3);
00527 hZShwVtx = new TH1D("hZShwVtx","",50,0,10);
00528
00529 if (event.simFlag == SimFlag::kMC){
00530
00531
00532 }
00533
00534 fvHistos.push_back(hRecoEnergyVsTime);
00535
00536 fvHistos.push_back(hRecoEnergyAll);
00537 fvHistos.push_back(hRecoEnergyNQ);
00538 fvHistos.push_back(hRecoEnergyPQ);
00539
00540
00541 fvHistos.push_back(hbatchPOT_1batchmode);
00542 fvHistos.push_back(hbatchPOT_2batchmode);
00543 fvHistos.push_back(hbatchPOT_3batchmode);
00544 fvHistos.push_back(hbatchPOT_4batchmode);
00545 fvHistos.push_back(hbatchPOT_5batchmode);
00546 fvHistos.push_back(hbatchPOT_6batchmode);
00547
00548 fvHistos.push_back(hCalibrateBatchPOT);
00549
00550 fvHistos.push_back(hBatchPOTSum_1batchmode);
00551 fvHistos.push_back(hBatchPOTSum_2batchmode);
00552 fvHistos.push_back(hBatchPOTSum_3batchmode);
00553 fvHistos.push_back(hBatchPOTSum_4batchmode);
00554 fvHistos.push_back(hBatchPOTSum_5batchmode);
00555 fvHistos.push_back(hBatchPOTSum_6batchmode);
00556
00557 fvHistos.push_back(hEventsPerTimePerPOTLT3GeV);
00558 fvHistos.push_back(hEventsPerTimePerPOTGT3GeV);
00559 fvHistos.push_back(hEventsPerTimePerPOTLT6GeV);
00560 fvHistos.push_back(hEventsPerTimePerPOTGT6GeV);
00561 fvHistos.push_back(hEventsPerTimePerPOTLT8GeV);
00562 fvHistos.push_back(hEventsPerTimePerPOTGT8GeV);
00563 fvHistos.push_back(hEventsPerTimePerPOTLT14GeV);
00564 fvHistos.push_back(hEventsPerTimePerPOTGT14GeV);
00565
00566
00567 fvHistos.push_back(hRecoEnergyAllmuon);
00568 fvHistos.push_back(hEventsPerTimePerPOTLT3GeVmuon);
00569 fvHistos.push_back(hEventsPerTimePerPOTGT3GeVmuon);
00570 fvHistos.push_back(hEventsPerTimePerPOTLT6GeVmuon);
00571 fvHistos.push_back(hEventsPerTimePerPOTGT6GeVmuon);
00572 fvHistos.push_back(hEventsPerTimePerPOTLT8GeVmuon);
00573 fvHistos.push_back(hEventsPerTimePerPOTGT8GeVmuon);
00574 fvHistos.push_back(hEventsPerTimePerPOTLT14GeVmuon);
00575 fvHistos.push_back(hEventsPerTimePerPOTGT14GeVmuon);
00576
00577 fvHistos.push_back(hSelectionBatchesAll);
00578 fvHistos.push_back(hSelectionBatches_BatchMode6);
00579 fvHistos.push_back(hSelectionBatches_BatchMode5);
00580 fvHistos.push_back(hSelectionBatches_BatchMode4);
00581 fvHistos.push_back(hSelectionBatches_BatchMode3);
00582 fvHistos.push_back(hSelectionBatches_BatchMode2);
00583 fvHistos.push_back(hSelectionBatches_BatchMode1);
00584 fvHistos.push_back(hVertexZvsEnergy);
00585
00586 fvHistos.push_back(hPotPerIntensity);
00587
00588 fvHistos.push_back(hPotPerIntensity_5BatchMode_Batch[0]);
00589 fvHistos.push_back(hPotPerIntensity_5BatchMode_Batch[1]);
00590 fvHistos.push_back(hPotPerIntensity_5BatchMode_Batch[2]);
00591 fvHistos.push_back(hPotPerIntensity_5BatchMode_Batch[3]);
00592 fvHistos.push_back(hPotPerIntensity_5BatchMode_Batch[4]);
00593
00594 fvHistos.push_back(hPotPerIntensity_6BatchMode_Batch[0]);
00595 fvHistos.push_back(hPotPerIntensity_6BatchMode_Batch[1]);
00596 fvHistos.push_back(hPotPerIntensity_6BatchMode_Batch[2]);
00597 fvHistos.push_back(hPotPerIntensity_6BatchMode_Batch[3]);
00598 fvHistos.push_back(hPotPerIntensity_6BatchMode_Batch[4]);
00599 fvHistos.push_back(hPotPerIntensity_6BatchMode_Batch[5]);
00600
00601 fvHistos.push_back(IntensityPerTime);
00602
00603 fvHistos.push_back(IntensityPerTime_5BatchMode_Batch[0]);
00604 fvHistos.push_back(IntensityPerTime_5BatchMode_Batch[1]);
00605 fvHistos.push_back(IntensityPerTime_5BatchMode_Batch[2]);
00606 fvHistos.push_back(IntensityPerTime_5BatchMode_Batch[3]);
00607 fvHistos.push_back(IntensityPerTime_5BatchMode_Batch[4]);
00608
00609 fvHistos.push_back(IntensityPerTime_6BatchMode_Batch[0]);
00610 fvHistos.push_back(IntensityPerTime_6BatchMode_Batch[1]);
00611 fvHistos.push_back(IntensityPerTime_6BatchMode_Batch[2]);
00612 fvHistos.push_back(IntensityPerTime_6BatchMode_Batch[3]);
00613 fvHistos.push_back(IntensityPerTime_6BatchMode_Batch[4]);
00614 fvHistos.push_back(IntensityPerTime_6BatchMode_Batch[5]);
00615
00616 fvHistos.push_back(hNumberofGoodSpillsPerIntensity);
00617
00618 fvHistos.push_back(hedgeActivityPH);
00619 fvHistos.push_back(hedgeActivityStrips);
00620 fvHistos.push_back(hminTimeSeparation);
00621 fvHistos.push_back(hoppEdgePH);
00622 fvHistos.push_back(hoppEdgeStrips);
00623 fvHistos.push_back(hcloseTimeDeltaZ);
00624 fvHistos.push_back(hnstripEvt);
00625 fvHistos.push_back(hstripsOverSqPlanes);
00626
00627 fvHistos.push_back(hydistribution);
00628 fvHistos.push_back(hplanesEvent);
00629 fvHistos.push_back(hPhEvent);
00630 fvHistos.push_back(hPhPerPlane);
00631 fvHistos.push_back(hPhPerStrip);
00632 fvHistos.push_back(hStripsPerPlane);
00633 fvHistos.push_back(hTracksInEvent);
00634 fvHistos.push_back(hShowersInEvent);
00635 fvHistos.push_back(hEventXVertex);
00636 fvHistos.push_back(hEventYVertex);
00637 fvHistos.push_back(hEventZVertex);
00638 fvHistos.push_back(hEventRVertex);
00639
00640 fvHistos.push_back(hTrkLength);
00641 fvHistos.push_back(hPlanesCrossedInTrk);
00642 fvHistos.push_back(hTrkCosTheta);
00643 fvHistos.push_back(hTrkdzds);
00644 fvHistos.push_back(hTrkPh);
00645 fvHistos.push_back(hTrkMomentumRange);
00646 fvHistos.push_back(hTrkPhPerPlane);
00647 fvHistos.push_back(hTrkqpSigma);
00648 fvHistos.push_back(hTrkPassedFit);
00649 fvHistos.push_back(hTrkChi2PerNdf);
00650 fvHistos.push_back(hTrkXVertex);
00651 fvHistos.push_back(hTrkYVertex);
00652 fvHistos.push_back(hTrkZVertex);
00653
00654 fvHistos.push_back(hPlanesCrossedInShower);
00655 fvHistos.push_back(hShowerEnergy);
00656 fvHistos.push_back(hXShwVtx);
00657 fvHistos.push_back(hYShwVtx);
00658 fvHistos.push_back(hZShwVtx);
00659
00660 fHistogramsCreated = true;
00661 }
00662
00664 void DataStabilityAnalysisClass::DeleteHistograms()
00665 {
00666 for (vector<TH1*>::iterator it = fvHistos.begin();
00667 it != fvHistos.end();
00668 ++it){
00669 this->Delete(*it);
00670 }
00671 }
00672
00674 void DataStabilityAnalysisClass::Delete(TH1* p)
00675 {
00676 if (p){delete p; p=0;}
00677 }
00678
00680 void DataStabilityAnalysisClass::FillTTree(NuEvent& event, BDSpillAccessor& bsAccess, BeamMonSpill* beamSpill, Bool_t passedSelectionCut)
00681 {
00682
00683
00684 if( event.goodBeamSntp && event.coilIsOk && ( event.isGoodDataQuality || !event.useDBForDataQuality ) ){
00685
00686
00687 ttimeRelToSpill = (event.timeEvtMin + event.crateT0);
00688 ttimeSec = event.timeSec;
00689 ttimeNanoSec = event.timeNanoSec;
00690 tnearestSpillSec = event.nearestSpillSec;
00691 tnearestSpillNanosec = event.nearestSpillNanosec;
00692 trun = event.run;
00693 tnshw = event.nshw;
00694 tsubRun = event.subRun;
00695 tsnarl = event.snarl;
00696 ttrkEn = event.trkEn;
00697 ttrkEnRange = event.trkEnRange;
00698 ttrkEnCurv = event.trkEnCurv;
00699 tshwEn = event.shwEn;
00700 tshwEnLinCCCor = event.shwEnCor;
00701 troID = event.roID;
00702
00703
00704 if(passedSelectionCut == true){
00705 tselectionevent = 1;
00706 }
00707 else{
00708 tselectionevent = 0;
00709 }
00710
00711
00712 VldTimeStamp timeStamp(event.timeSec,event.timeNanoSec);
00713
00714 if( event.snarl != prevSnarl ) {
00715
00716 DSTspillcounter++;
00717
00718 beamSpill = const_cast<BeamMonSpill*>( bsAccess.LoadSpill(timeStamp) );
00719 VldTimeStamp vldBeamTime = beamSpill->SpillTime();
00720
00721 beamTimeSec = (ULong64_t)vldBeamTime.GetSec();
00722 beamTimeNSec = (ULong64_t)vldBeamTime.GetNanoSec();
00723
00724 fTor101 = beamSpill->fTor101;
00725
00726 fTr101d = beamSpill->fTr101d;
00727 fTortgt = beamSpill->fTortgt;
00728 fTrtgtd = beamSpill->fTrtgtd;
00729
00730 DSTpotCounter += fTrtgtd;
00731
00732 fHornCur = beamSpill->fHornCur;
00733
00734
00735
00736 memcpy( fTargBpmX, beamSpill->fTargBpmX, 6*sizeof(float) );
00737 memcpy( fTargBpmY, beamSpill->fTargBpmY, 6*sizeof(float) );
00738 memcpy( fBpmInt, beamSpill->fBpmInt, 6*sizeof(float) );
00739
00740 fTargProfX = beamSpill->fTargProfX;
00741 fTargProfY = beamSpill->fTargProfY;
00742 fProfWidX = beamSpill->fProfWidX;
00743
00744 fProfWidY = beamSpill->fProfWidY;
00745 fHadInt = beamSpill->fHadInt;
00746
00747 fMuInt1 = beamSpill->fMuInt1;
00748 fMuInt2 = beamSpill->fMuInt2;
00749 fMuInt3 = beamSpill->fMuInt3;
00750 totBpm = fBpmInt[0] + fBpmInt[1] + fBpmInt[2] + fBpmInt[3] + fBpmInt[4] + fBpmInt[5];
00751
00752 numBatches = 0;
00753
00754 for( int batch = 0; batch < 6; batch++ ){
00755
00756 if( fBpmInt[batch] > 0 ){
00757 numBatches++;
00758 }
00759 else{
00760
00761 }
00762 }
00763
00764 hCalibrateBatchPOT->Fill(fTrtgtd, totBpm);
00765
00766 batchPot[0] = 0;
00767 batchPot[1] = 0;
00768 batchPot[2] = 0;
00769 batchPot[3] = 0;
00770 batchPot[4] = 0;
00771 batchPot[5] = 0;
00772
00773 for(int batchnum = 0; batchnum < 6; ++batchnum){
00774
00775 batchPot[batchnum] = fBpmInt[batchnum] * fTrtgtd / totBpm;
00776
00777 if(numBatches == 1){
00778 hbatchPOT_1batchmode->Fill(batchnum, batchPot[batchnum]);
00779 hBatchPOTSum_1batchmode->Fill(batchnum, batchPot[batchnum]);
00780 }
00781 else if(numBatches == 2){
00782 hbatchPOT_2batchmode->Fill(batchnum, batchPot[batchnum]);
00783 hBatchPOTSum_2batchmode->Fill(batchnum, batchPot[batchnum]);
00784 }
00785 else if(numBatches == 3){
00786 hbatchPOT_3batchmode->Fill(batchnum, batchPot[batchnum]);
00787 hBatchPOTSum_3batchmode->Fill(batchnum, batchPot[batchnum]);
00788 }
00789 else if(numBatches == 4){
00790 hbatchPOT_4batchmode->Fill(batchnum, batchPot[batchnum]);
00791 hBatchPOTSum_4batchmode->Fill(batchnum, batchPot[batchnum]);
00792 }
00793 else if(numBatches == 5){
00794 hbatchPOT_5batchmode->Fill(batchnum, batchPot[batchnum]);
00795 hBatchPOTSum_5batchmode->Fill(batchnum, batchPot[batchnum]);
00796 }
00797 else if(numBatches == 6){
00798 hbatchPOT_6batchmode->Fill(batchnum, batchPot[batchnum]);
00799 hBatchPOTSum_6batchmode->Fill(batchnum, batchPot[batchnum]);
00800 if(batchnum == 0) {batchBin1summer += batchPot[batchnum];}
00801 }
00802 else{
00803 assert(false);
00804
00805 }
00806 }
00807
00808 cout << "snarl " << event.snarl << "Run number = " << event.run << ", POT 6 batch Bin1 = " << batchBin1summer << endl;
00809
00810
00811 this->IntensityProfile(event, batchPot, fTrtgtd, numBatches);
00812
00813
00814
00815
00816
00817
00818
00819
00820
00821
00822 }
00823
00824
00825 if(numBatches == 6) hSelectionBatches_BatchMode6->Fill((event.timeEvtMin + event.crateT0), event.rw);
00826 else if (numBatches == 5) hSelectionBatches_BatchMode5->Fill((event.timeEvtMin + event.crateT0), event.rw);
00827 else if (numBatches == 4) hSelectionBatches_BatchMode4->Fill((event.timeEvtMin + event.crateT0), event.rw);
00828 else if (numBatches == 3) hSelectionBatches_BatchMode3->Fill((event.timeEvtMin + event.crateT0), event.rw);
00829 else if (numBatches == 2) hSelectionBatches_BatchMode2->Fill((event.timeEvtMin + event.crateT0), event.rw);
00830 else if (numBatches == 1) hSelectionBatches_BatchMode1->Fill((event.timeEvtMin + event.crateT0), event.rw);
00831 else{
00832
00833 }
00834
00835
00836 whichBatch=-1;
00837
00838
00839 if(ttimeRelToSpill > 1.28e-6 && ttimeRelToSpill <= 2.9e-6 ) whichBatch = 0;
00840 else if(ttimeRelToSpill > 2.9e-6 && ttimeRelToSpill <= 4.49e-6 ) whichBatch = 1;
00841 else if(ttimeRelToSpill > 4.49e-6 && ttimeRelToSpill <= 6.06e-6 ) whichBatch = 2;
00842 else if(ttimeRelToSpill > 6.06e-6 && ttimeRelToSpill <= 7.65e-6 ) whichBatch = 3;
00843 else if(ttimeRelToSpill > 7.65e-6 && ttimeRelToSpill <= 9.23e-6 ) whichBatch = 4;
00844 else if(ttimeRelToSpill > 9.23e-6 && ttimeRelToSpill <= 10.82e-6) whichBatch = 5;
00845
00846 timeTree->Fill();
00847 }
00848
00849 prevSnarl = event.snarl;
00850
00851 }
00853 void DataStabilityAnalysisClass::FillPlots(NuEvent& event)
00854 {
00855 this->PrepareEventForPlotting(event);
00856
00857 Double_t energy = 0.0;
00858
00859
00860
00861 if(fCutter->AnaVersion() == NuCuts::kCC1070Std){
00862 energy = event.energy;
00863 }
00864 if(fCutter->AnaVersion() == NuCuts::kNCPRLNC){
00865 energy = event.shwEnLinNCCor;
00866 }
00867
00868
00869 hRecoEnergyVsTime->Fill(FixedTime, energy);
00870
00871 hRecoEnergyAll->Fill( energy, 1);
00872
00873 if ( -1 == event.charge ){ hRecoEnergyNQ->Fill(energy, event.rw); }
00874 if ( 1 == event.charge ){ hRecoEnergyPQ->Fill(energy, event.rw); }
00875
00876
00877 if ( energy < 8.0 ){ hEventsPerTimePerPOTLT8GeV->Fill( FixedTime ); }
00878 if ( energy >= 8.0 ){ hEventsPerTimePerPOTGT8GeV->Fill( FixedTime ); }
00879 if ( energy < 6.0 ){ hEventsPerTimePerPOTLT6GeV->Fill( FixedTime ); }
00880 if ( energy >= 6.0 ){ hEventsPerTimePerPOTGT6GeV->Fill( FixedTime ); }
00881 if ( energy < 3.0 ){ hEventsPerTimePerPOTLT3GeV->Fill( FixedTime ); }
00882 if ( energy >= 3.0 ){ hEventsPerTimePerPOTGT3GeV->Fill( FixedTime ); }
00883 if ( energy < 14.0){ hEventsPerTimePerPOTLT14GeV->Fill( FixedTime); }
00884 if ( energy >= 14.0){ hEventsPerTimePerPOTGT14GeV->Fill( FixedTime); }
00885
00886
00887 if(fCutter->AnaVersion()==29){
00888
00889 if ( event.trkEn < 8.0){ hEventsPerTimePerPOTLT8GeVmuon->Fill( FixedTime ); }
00890 if ( event.trkEn >= 8.0){ hEventsPerTimePerPOTGT8GeVmuon->Fill( FixedTime ); }
00891 if ( event.trkEn < 6.0){ hEventsPerTimePerPOTLT6GeVmuon->Fill( FixedTime ); }
00892 if ( event.trkEn >= 6.0){ hEventsPerTimePerPOTGT6GeVmuon->Fill( FixedTime ); }
00893 if ( event.trkEn < 3.0){ hEventsPerTimePerPOTLT3GeVmuon->Fill( FixedTime ); }
00894 if ( event.trkEn >= 3.0){ hEventsPerTimePerPOTGT3GeVmuon->Fill( FixedTime ); }
00895 if ( event.trkEn < 14.0){hEventsPerTimePerPOTLT14GeVmuon->Fill( FixedTime); }
00896 if ( event.trkEn >= 14.0){hEventsPerTimePerPOTGT14GeVmuon->Fill( FixedTime); }
00897
00898 hRecoEnergyAllmuon->Fill(event.trkEn, event.rw);
00899 }
00900
00901
00902
00903 hedgeActivityPH->Fill( event.edgeActivityPH);
00904 hedgeActivityStrips->Fill(event.edgeActivityStrips);
00905 hminTimeSeparation->Fill( event.minTimeSeparation*1e9);
00906 hoppEdgePH->Fill( event.oppEdgePH);
00907 hoppEdgeStrips->Fill( event.oppEdgeStrips);
00908 hcloseTimeDeltaZ->Fill( event.closeTimeDeltaZ);
00909 hnstripEvt->Fill( event.nstripEvt);
00910
00911 Double_t stripsPerSqPlanes = 0;
00912 stripsPerSqPlanes = ((Double_t)event.nstripEvt)/(event.planeEvtN*event.planeEvtN);
00913 hstripsOverSqPlanes->Fill( stripsPerSqPlanes);
00914
00915
00916 hydistribution->Fill(event.y, event.rw);
00917 hplanesEvent->Fill( event.planeEvtN);
00918 hPhEvent->Fill( event.rawPhEvt);
00919 hPhPerPlane->Fill( ((Double_t)event.rawPhEvt/(Double_t)event.planeEvtN));
00920 hPhPerStrip->Fill( (Double_t)(event.rawPhEvt/(Double_t)event.nstripEvt));
00921 hStripsPerPlane->Fill((Double_t)(event.nstripEvt/event.planeEvtN));
00922 hTracksInEvent->Fill( event.ntrk);
00923 hShowersInEvent->Fill(event.nshw);
00924 hEventXVertex->Fill( event.xEvtVtx);
00925 hEventYVertex->Fill( event.yEvtVtx);
00926 hEventZVertex->Fill( event.zEvtVtx);
00927
00928 hSelectionBatchesAll->Fill((event.timeEvtMin + event.crateT0), event.rw);
00929 hVertexZvsEnergy ->Fill(energy, event.zEvtVtx, event.rw);
00930
00931
00932 Double_t RVertex = TMath::Sqrt((event.xEvtVtx * event.xEvtVtx) +
00933 (event.yEvtVtx * event.yEvtVtx));
00934 hEventRVertex->Fill(RVertex);
00935
00936
00937 hTrkLength->Fill( event.trkLength);
00938 hPlanesCrossedInTrk->Fill(event.trknplane);
00939 hTrkCosTheta->Fill( event.dirCosNu);
00940
00941 hTrkdzds->Fill( event.trkvtxdcosz);
00942 hTrkPh->Fill( event.trkphsigcor);
00943 hTrkMomentumRange->Fill( event.trkMomentumRange);
00944 hTrkPhPerPlane->Fill( event.knn10TrkMeanPh);
00945 hTrkqpSigma->Fill( event.qp_sigqp);
00946 hTrkPassedFit->Fill( event.trkfitpass);
00947 hTrkChi2PerNdf->Fill( event.chi2PerNdof);
00948 hTrkXVertex->Fill( event.xTrkVtx);
00949 hTrkYVertex->Fill( event.yTrkVtx);
00950 hTrkZVertex->Fill( event.zTrkVtx);
00951
00952
00953 hPlanesCrossedInShower->Fill(event.nplaneShw);
00954 hShowerEnergy->Fill(event.shwEn);
00955 hXShwVtx->Fill(event.xShwVtx);
00956 hYShwVtx->Fill(event.yShwVtx);
00957 hZShwVtx->Fill(event.zShwVtx);
00958
00959 }
00960
00962 void DataStabilityAnalysisClass::WritePlots()
00963 {
00964 for (vector<TH1*>::iterator it = fvHistos.begin();
00965 it != fvHistos.end();
00966 ++ it){
00967 (*it)->Write();
00968 }
00969
00970
00971 timeTree->Write("timeTree");
00972
00973
00974 std::cout << " " << std::endl;
00975 std::cout << "Spill POT = " << DSTpotCounter << std::endl;
00976 std::cout << "Spill counter = " << DSTspillcounter << std::endl;
00977 std::cout << " " << std::endl;
00978 }
00980 void DataStabilityAnalysisClass::IntensityProfile(NuEvent& nu, Double_t* batchPot, Double_t intensity, Int_t numBatches)
00981 {
00982
00983
00984
00985
00986
00987
00988 hPotPerIntensity ->Fill( intensity, intensity );
00989 hNumberofGoodSpillsPerIntensity->Fill( intensity );
00990 IntensityPerTime ->Fill( FixedTime, intensity, 1);
00991
00992
00993 if(numBatches == 5){
00994 for(int i = 0; i < 5; ++i){
00995 hPotPerIntensity_5BatchMode_Batch[i]->Fill(batchPot[i]);
00996 IntensityPerTime_5BatchMode_Batch[i]->Fill(FixedTime, batchPot[i], 1);
00997 }
00998 }
00999 else if(numBatches == 6){
01000 for(int i = 0; i < 6; ++i){
01001 hPotPerIntensity_6BatchMode_Batch[i]->Fill(batchPot[i]);
01002 IntensityPerTime_6BatchMode_Batch[i]->Fill(FixedTime, batchPot[i], 1);
01003 }
01004 }
01005 else{
01006
01007 }
01008
01009 }
01011 void DataStabilityAnalysisClass::ObtainPotPerRun()
01012 {
01013 std::cout << "\n----------------------Obtaining Pot Per Run-----------------" << std::endl;
01014
01015 string inputFileName = this->GetInputFileName();
01016
01017 std::cout << "The InputFilename is: " << this->GetInputFileName() << endl;
01018
01019 std::cout << "Creating the TChain" << std::endl;
01020
01021
01022 TChain tChain("s");
01023 tChain.Add(inputFileName.c_str());
01024 TObjArray* fileList = tChain.GetListOfFiles();
01025
01026 TH1D* hRunPot = new TH1D("hRunPot", "hRunPot; Run; Pot", 20000, 7000, 27000);
01027 TH1D* hGoodRunSpills = new TH1D("hGoodRunSpills", "hGoodRunSpills; Run; # Good Spills", 20000, 7000, 27000);
01028 TH1D* hTotalRunSpills = new TH1D("hTotalRunSpills", "hTotalRunSpills; Run; # Total Spills", 20000, 7000, 27000);
01029 TH1D* hPotPerTime = new TH1D("hPotPerTime", "hPotPerTime", 3652, StartTime, EndTime);
01030
01031 TH1D* hTotalSnarlCount = new TH1D("TotalSnarlCount","TotalSnarlCount",2,0,2);
01032
01033 TH1D* hTotalBwidx = new TH1D("hTotalBwidx", "hTotalBwidx; beam width x; # of spills", 1000, -2e-3, 10e-3);
01034 TH1D* hTotalBwidy = new TH1D("hTotalBwidy", "hTotalBwidy; beam width y; # of spills", 1000, -2e-3, 10e-3);
01035 TH1D* hTotalBposx = new TH1D("hTotalBposx", "hTotalBposx; beam pos x; # of spills", 5000, -10e-3, 10e-3);
01036 TH1D* hTotalBposy = new TH1D("hTotalBposy", "hTotalBposy; beam pos y; # of spills", 5000, -10e-3, 10e-3);
01037
01038 TH2D* TimeVswidx = new TH2D("TimeVswidx","TimeVswidx; beam width x; time (s)", 500, -0.002, 0.002, 3652, StartTime, EndTime);
01039 TH2D* TimeVswidy = new TH2D("TimeVswidy","TimeVswidy; beam width y; time (s)", 500, -0.002, 0.002, 3652, StartTime, EndTime);
01040 TH2D* TimeVsposx = new TH2D("TimeVsposx","TimeVsposx; beam pos x; time (s)", 500, 0.0005, 0.0013, 3652, StartTime, EndTime);
01041 TH2D* TimeVsposy = new TH2D("TimeVsposy","TimeVsposy; beam pos y; time (s)", 500, 0.0014, 0.0022, 3652, StartTime, EndTime);
01042
01043 NuEvent *event(0);
01044
01045
01046 Int_t TotalSnarlCount = -999;
01047
01048
01049 TTimeStamp* s = new TTimeStamp();
01050
01051 std::cout << "There are " << fileList->GetEntries() << " root files to loop over." << std::endl;
01052
01053 for ( Int_t counter=0; counter<fileList->GetEntries(); ++counter ){
01054
01055 TChainElement* chainElement = dynamic_cast<TChainElement*> (fileList->At(counter));
01056
01057 if ( chainElement ){
01058
01059 string fileName = chainElement->GetTitle();
01060 TFile f(fileName.c_str(),"READ");
01061
01062 cout << "Filename = " << fileName.c_str() << endl;
01063
01064 TTree *tree= (TTree*)f.Get("s");
01065 tree->SetBranchAddress("s",&event);
01066 tree->GetEntry(0);
01067
01068
01069 TH1D* hTotalPot = (TH1D*)f.Get("hTotalPot"); assert(hTotalPot);
01070 TH1D* hSnarlCount = (TH1D*)f.Get("hSnarlCount"); assert(hSnarlCount);
01071 TH1D* hBwidx = (TH1D*)f.Get("hBwidx"); assert(hBwidx);
01072 TH1D* hBwidy = (TH1D*)f.Get("hBwidy"); assert(hBwidy);
01073 TH1D* hBposx = (TH1D*)f.Get("hBposx"); assert(hBposx);
01074 TH1D* hBposy = (TH1D*)f.Get("hBposy"); assert(hBposy);
01075
01076 std::cout << "The Pot for file " << fileName.c_str() << " is " << hTotalPot->Integral() << std::endl;
01077 hPotPerTime->Fill( event->timeSeconds, hTotalPot->Integral() );
01078
01079 TotalSnarlCount += hSnarlCount->Integral();
01080
01081
01082 s->SetSec(event->timeSeconds);
01083
01084 cout << std::setprecision(10) << "the time = " << event->timeSeconds << endl;
01085
01086
01087 s->Print();
01088
01089 hRunPot->Fill( event->run, hTotalPot->Integral() );
01090
01091
01092 hGoodRunSpills->Fill( event->run, hBwidx->GetEntries() );
01093 hTotalRunSpills->Fill( event->run, hSnarlCount->Integral() );
01094
01095 hTotalBwidx->Add(hBwidx);
01096 hTotalBwidy->Add(hBwidy);
01097 hTotalBposx->Add(hBposx);
01098 hTotalBposy->Add(hBposy);
01099
01100 for(Int_t bin=1; bin<=hBwidx->GetNbinsX(); ++bin){
01101 TimeVswidx->Fill(hBwidx->GetXaxis()->GetBinCenter(bin), event->timeSeconds, hBwidx->GetBinContent(bin));
01102 }
01103 for(Int_t bin=1; bin<=hBwidy->GetNbinsX(); ++bin){
01104 TimeVswidy->Fill(hBwidy->GetXaxis()->GetBinCenter(bin), event->timeSeconds, hBwidy->GetBinContent(bin));
01105 }
01106 for(Int_t bin=1; bin<=hBposx->GetNbinsX(); ++bin){
01107 TimeVsposx->Fill(hBposx->GetXaxis()->GetBinCenter(bin), event->timeSeconds, hBposx->GetBinContent(bin));
01108 }
01109 for(Int_t bin=1; bin<=hBposy->GetNbinsX(); ++bin){
01110 TimeVsposy->Fill(hBposy->GetXaxis()->GetBinCenter(bin), event->timeSeconds, hBposy->GetBinContent(bin));
01111 }
01112
01113 f.Close();
01114 }
01115 }
01116
01117 std::cout << " " << std::endl;
01118 std::cout << "TotalSnarlCount for this set of files = " << TotalSnarlCount << std::endl;
01119 std::cout << " " << std::endl;
01120
01121 hTotalSnarlCount->Fill(1, TotalSnarlCount);
01122
01123 fOutFile->cd();
01124
01125 hRunPot->Write();
01126 hTotalSnarlCount->Write();
01127 hTotalRunSpills->Write();
01128 hGoodRunSpills->Write();
01129 hPotPerTime->Write();
01130 hTotalBwidx->Write();
01131 hTotalBwidy->Write();
01132 hTotalBposx->Write();
01133 hTotalBposy->Write();
01134 TimeVswidx->Write();
01135 TimeVswidy->Write();
01136 TimeVsposx->Write();
01137 TimeVsposy->Write();
01138
01139 std::cout << " " << std::endl;
01140
01141 }
01142
01143
01144
01145