BOSS 7.1.2
BESIII Offline Software System
Loading...
Searching...
No Matches
MucCalibMgr.cxx
Go to the documentation of this file.
1//------------------------------------------------------------------------------|
2// [File ]: MucCalibMgr.cxx |
3// [Brief ]: Manager of MUC calibration algrithom |
4// [Author]: Xie Yuguang, <[email protected]> |
5// [Date ]: Mar 28, 2006 |
6// [Log ]: See ChangLog |
7//------------------------------------------------------------------------------|
8#include <cstdio>
9#include <iostream>
10#include <fstream>
11
12#include "CLHEP/Vector/LorentzVector.h"
13#include "CLHEP/Vector/ThreeVector.h"
14#include "CLHEP/Vector/LorentzVector.h"
15#include "CLHEP/Vector/TwoVector.h"
16#include "CLHEP/Geometry/Point3D.h"
17using CLHEP::Hep3Vector;
18using CLHEP::Hep2Vector;
19using CLHEP::HepLorentzVector;
20
21#include "GaudiKernel/IHistogramSvc.h"
22#include "GaudiKernel/MsgStream.h"
23#include "GaudiKernel/AlgFactory.h"
24#include "GaudiKernel/ISvcLocator.h"
25#include "GaudiKernel/SmartDataPtr.h"
26#include "GaudiKernel/SmartDataLocator.h"
27#include "GaudiKernel/IDataProviderSvc.h"
28#include "GaudiKernel/Bootstrap.h"
29#include "GaudiKernel/PropertyMgr.h"
30#include "GaudiKernel/INTupleSvc.h"
31#include "GaudiKernel/NTuple.h"
32#include "GaudiKernel/Algorithm.h"
33
34#include "Identifier/Identifier.h"
35#include "Identifier/MucID.h"
36#include "EventModel/EventHeader.h"
37#include "EventModel/EventModel.h"
38#include "EvTimeEvent/RecEsTime.h"
39
40#include "McTruth/McKine.h"
41#include "McTruth/McParticle.h"
42#include "McTruth/MucMcHit.h"
43
44#include "MdcRecEvent/RecMdcTrack.h"
45#include "EmcRecEventModel/RecEmcEventModel.h"
46#include "TofRecEvent/RecTofTrack.h"
47#include "DstEvent/TofHitStatus.h"
48
49//#include "ReconEvent/ReconEvent.h"
50#include "MucRawEvent/MucDigi.h"
51#include "MucRecEvent/MucRecHit.h"
52#include "MucRecEvent/RecMucTrack.h"
53#include "MucRecEvent/MucRecHitContainer.h"
54
55#include "MucCalibAlg/MucStructConst.h"
56#include "MucCalibAlg/MucCalibMgr.h"
57#include "MucCalibAlg/MucMark.h"
58#include "MucCalibAlg/MucBoxCal.h"
59#include "MucCalibAlg/MucStripCal.h"
60
61using namespace std;
62using namespace Event;
63
64//------------------------------------------- Constructor -----------------------------------------
65MucCalibMgr::MucCalibMgr( std::vector<double> jobInfo, std::vector<int> configInfo, std::string outputFile )
66{
67 MsgStream log(msgSvc, "MucCalibMgr");
68 log << MSG::INFO << "-------Initializing------- " << endreq;
69
70 // init Gaudi service
71 log << MSG::INFO << "Initialize Gaudi service" << endreq;
72 Gaudi::svcLocator() -> service("MessageSvc", msgSvc);
73 ntupleSvc = NULL;
74 eventSvc = NULL;
75
76 m_jobStart = clock();
77
78 // init parameters
79 log << MSG::INFO << "Initialize parameters" << endreq;
80
81 m_vs = jobInfo[0]; m_hv = jobInfo[1]; m_th = jobInfo[2];
82
83 if( jobInfo[3] <= 0 ) m_er = TRIGGER_RATE; // 4000Hz
84 else m_er = jobInfo[3];
85
86 if( jobInfo[4] <= 0 || jobInfo[4] >5 ) m_tag = 0;
87 else m_tag = jobInfo[4];
88
89 if( configInfo[0] < 0 || configInfo[0] > 2 )
90 m_recMode = 0;
91 else
92 m_recMode = configInfo[0];
93
94 m_usePad = configInfo[1];
95
96 if( configInfo[2] < 0 || configInfo[2] > STRIP_INBOX_MAX )
97 m_effWindow = EFF_WINDOW; // 4 strip-width
98 else m_effWindow = configInfo[2];
99
100 if( configInfo[3]<0 || configInfo[3] >4 )
101 m_clusterMode = DEFAULT_BUILD_MODE;
102 else m_clusterMode = configInfo[3];
103
104 m_clusterSave = configInfo[4];
105 m_checkEvent = configInfo[5];
106 m_dimuSelect = configInfo[6];
107 m_dimuOnly = configInfo[7];
108
109 m_outputFile = outputFile;
110
111 m_fdata = NULL;
112 if( m_clusterMode ==1 && m_clusterSave == 1 )
113 m_fdata = new ofstream("MucCluster.dat", ios::out);
114
115 m_fStartRun = m_fEndRun = 0;
116 m_currentRun = m_currentEvent = 0;
117 m_eventTag = 0;
118
119 for( int i=0; i<PART_MAX; i++ )
120 for( int j=0; j<SEGMENT_MAX; j++ )
121 for( int k=0; k<LAYER_MAX; k++ )
122 for( int n=0; n<STRIP_INBOX_MAX; n++ )
123 {
124 m_record[i][j][k][n][0] = 0;
125 m_record[i][j][k][n][1] = 0;
126 m_record[i][j][k][n][2] = 0;
127 m_record[i][j][k][n][3] = 0;
128 m_record[i][j][k][n][4] = 0;
129 }
130
131 for( int i=0; i<LAYER_MAX; i++ )
132 {
133 m_layerResults[0][i] = 0;
134 m_layerResults[1][i] = 0;
135 m_layerResults[2][i] = 0;
136 m_layerResults[3][i] = 0;
137 m_layerResults[4][i] = 0;
138 }
139
140 for( int i=0; i<BOX_MAX; i++ )
141 {
142 m_boxResults[0][i] = 0;
143 m_boxResults[1][i] = 0;
144 m_boxResults[2][i] = 0;
145 m_boxResults[3][i] = 0;
146 m_boxResults[4][i] = 0;
147 }
148
149 for( int i=0; i<STRIP_MAX; i++ )
150 {
151 m_stripResults[0][i] = 0;
152 m_stripResults[1][i] = 0;
153 m_stripResults[2][i] = 0;
154 m_stripResults[3][i] = 0;
155 // strip no cluster
156 }
157
158 m_fTotalDAQTime = 0;
159
160 m_ptrMucMark = new MucMark(0,0,0,0);
161 m_ptrIdTr = new MucIdTransform();
162
163 m_digiCol.resize(0);
164 m_calHitCol.resize(0);
165 m_expHitCol.resize(0);
166 m_effHitCol.resize(0);
167 m_nosHitCol.resize(0);
168 m_clusterCol.resize(0);
169
170 for( int i=0; i<PART_MAX; i++ )
171 for( int j=0; j<SEGMENT_MAX; j++ ) {
172 m_segDigiCol[i][j].resize(0);
173 }
174
175 // init Gaudi Ntuple
176 InitNtuple();
177
178 // init ROOT histograms
179 InitHisto();
180
181 // init Tree
182 InitConstTree();
183
184 // init area of units
185 InitArea();
186 log << MSG::INFO << "Initialize done!" << endreq;
187}
188
189//----------------------------------------------- Destructor --------------------------------------
190MucCalibMgr::~MucCalibMgr()
191{
192 ClearOnlineHisto();
193
194 ClearHistoLV2();
195 ClearHistoLV1();
196 ClearHistoLV0();
197 if( m_usePad != 0 ) Clear2DEffHisto();
198
199 ClearClusterHisto();
200 ClearResHisto();
201 ClearConstTree();
202
203 if(m_ptrMucMark) delete m_ptrMucMark;
204 if(m_ptrIdTr) delete m_ptrIdTr;
205
206 m_digiCol.clear();
207 m_calHitCol.clear();
208 m_effHitCol.clear();
209 m_nosHitCol.clear();
210 m_clusterCol.clear();
211
212 for( int i=0; i<PART_MAX; i++ )
213 for( int j=0; j<SEGMENT_MAX; j++ ) {
214 m_segDigiCol[i][j].clear();
215 }
216
217 if(m_record) delete []m_record;
218 if(m_layerResults) delete []m_layerResults;
219 if(m_boxResults) delete []m_boxResults;
220 if(m_stripResults) delete []m_stripResults;
221}
222
223// Clear online histograms
224StatusCode MucCalibMgr::ClearOnlineHisto()
225{
226 if(m_hHitMapBarrel_Lay) delete []m_hHitMapBarrel_Lay;
227 if(m_hHitMapEndcap_Lay) delete []m_hHitMapEndcap_Lay;
228 if(m_hHitMapBarrel_Seg) delete []m_hHitMapBarrel_Seg;
229 if(m_hHitMapEndcap_Seg) delete []m_hHitMapEndcap_Seg;
230
231 if(m_hHitVsEvent) delete m_hHitVsEvent;
232 if(m_hTrackDistance) delete m_hTrackDistance;
233 if(m_hTrackPosPhiDiff) delete m_hTrackPosPhiDiff;
234 if(m_hTrackPosThetaDiff)delete m_hTrackPosThetaDiff;
235 if(m_hTrackMomPhiDiff) delete m_hTrackMomPhiDiff;
236 if(m_hTrackMomThetaDiff)delete m_hTrackMomThetaDiff;
237 if(m_hDimuTracksPosDiff)delete m_hDimuTracksPosDiff;
238 if(m_hDimuTracksMomDiff)delete m_hDimuTracksMomDiff;
239 if(m_hPhiCosTheta) delete m_hPhiCosTheta;
240
241 return StatusCode::SUCCESS;
242}
243
244// Clear level 0 histograms
245StatusCode MucCalibMgr::ClearHistoLV0()
246{
247 if(m_hBrLayerFire) delete m_hBrLayerFire;
248 if(m_hEcLayerFire) delete m_hEcLayerFire;
249 if(m_hLayerFire) delete m_hLayerFire;
250 if(m_hLayerExpHit) delete m_hLayerExpHit;
251 if(m_hLayerExpHit) delete m_hLayerEffHit;
252 if(m_hLayerNosHit) delete m_hLayerNosHit;
253 if(m_hLayerEff) delete m_hLayerEff;
254 if(m_hLayerNosRatio)delete m_hLayerNosRatio;
255 if(m_hLayerArea) delete m_hLayerArea;
256 if(m_hLayerNos) delete m_hLayerNos;
257 if(m_hLayerCnt) delete m_hLayerCnt;
258
259 return StatusCode::SUCCESS;
260}
261
262// Clear level 1 histograms
263StatusCode MucCalibMgr::ClearHistoLV1()
264{
265 if(m_hBoxFire) delete m_hBoxFire;
266 if(m_hBoxExpHit) delete m_hBoxExpHit;
267 if(m_hBoxEffHit) delete m_hBoxEffHit;
268 if(m_hBoxNosHit) delete m_hBoxNosHit;
269 if(m_hBoxEff) delete m_hBoxEff;
270 if(m_hBoxNosRatio)delete m_hBoxNosRatio;
271 if(m_hBoxArea) delete m_hBoxArea;
272 if(m_hBoxNos) delete m_hBoxNos;
273 if(m_hBoxCnt) delete m_hBoxCnt;
274
275 return StatusCode::SUCCESS;
276}
277
278// Clear level 2 histograms
279StatusCode MucCalibMgr::ClearHistoLV2()
280{
281 for( int i=0; i< BOX_MAX; i++ )
282 {
283 if(m_hStripFireMap[i]) delete m_hStripFireMap[i];
284 if(m_hStripExpHitMap[i]) delete m_hStripExpHitMap[i];
285 if(m_hStripEffHitMap[i]) delete m_hStripEffHitMap[i];
286 if(m_hStripNosHitMap[i]) delete m_hStripNosHitMap[i];
287 if(m_hStripEffMap[i]) delete m_hStripEffMap[i];
288 if(m_hStripNosRatioMap[i])delete m_hStripNosRatioMap[i];
289 }
290
291 if(m_hStripFire) delete m_hStripFire;
292 if(m_hStripExpHit) delete m_hStripExpHit;
293 if(m_hStripEffHit) delete m_hStripEffHit;
294 if(m_hStripNosHit) delete m_hStripNosHit;
295 if(m_hStripEff) delete m_hStripEff;
296 if(m_hStripNosRatio)delete m_hStripNosRatio;
297 if(m_hStripArea) delete m_hStripArea;
298 if(m_hStripNos) delete m_hStripNos;
299 if(m_hStripCnt) delete m_hStripCnt;
300
301 return StatusCode::SUCCESS;
302}
303
304// Clear 2D eff histogrames
305StatusCode MucCalibMgr::Clear2DEffHisto()
306{
307
308 for( int i=0; i<PART_MAX; i++ )
309 for( int j=0; j<((i==BRID)?B_SEG_NUM:E_SEG_NUM); j++ )
310 for( int k=0; k<((i==BRID)?B_LAY_NUM:E_LAY_NUM); k++ )
311 {
312 if(m_h2DExpMap[i][j][k]) delete m_h2DExpMap[i][j][k];
313 if(m_h2DHitMap[i][j][k]) delete m_h2DHitMap[i][j][k];
314 if(m_h2DEffMap[i][j][k]) delete m_h2DEffMap[i][j][k];
315 }
316
317 return StatusCode::SUCCESS;
318}
319
320// Clear cluster histograms
321StatusCode MucCalibMgr::ClearClusterHisto()
322{
323 for( int i=0; i<LAYER_MAX; i++ ) {
324 if(m_hLayerCluster[i]) delete m_hLayerCluster[i];
325 }
326
327 for( int i=0; i<BOX_MAX; i++ ) {
328 if(m_hBoxCluster[i]) delete m_hBoxCluster[i];
329 }
330
331 if(m_hLayerClusterCmp) delete m_hLayerClusterCmp;
332 if(m_hBoxClusterCmp) delete m_hBoxClusterCmp;
333
334 return StatusCode::SUCCESS;
335}
336
337// Clear resolution histograms
338StatusCode MucCalibMgr::ClearResHisto()
339{
340 for( int i=0; i<B_LAY_NUM; i++ ) {
341 if(m_hBarrelResDist[i]) delete m_hBarrelResDist[i];
342 }
343 for( int i=0; i<E_LAY_NUM; i++ ) {
344 if(m_hEndcapResDist[i]) delete m_hEndcapResDist[i];
345 }
346
347 if(m_hBarrelResComp[0]) delete m_hBarrelResComp[0];
348 if(m_hBarrelResComp[1]) delete m_hBarrelResComp[1];
349 if(m_hEndcapResComp[0]) delete m_hEndcapResComp[0];
350 if(m_hEndcapResComp[1]) delete m_hEndcapResComp[1];
351
352 return StatusCode::SUCCESS;
353}
354
355// Clear Tree
356StatusCode MucCalibMgr::ClearConstTree()
357{
358
359 if(m_tJobLog) delete m_tJobLog;
360 if(m_tStatLog) delete m_tStatLog;
361 if(m_tLayConst) delete m_tLayConst;
362 if(m_tBoxConst) delete m_tBoxConst;
363 if(m_tStrConst) delete m_tStrConst;
364
365 return StatusCode::SUCCESS;
366}
367
368//-------------------------------------------------------------------------------------------------
369//-------------------------------- Member functions for initializing ------------------------------
370//-------------------------------------------------------------------------------------------------
371
372
373// init Ntuples
374StatusCode MucCalibMgr::InitNtuple()
375{
376 MsgStream log(msgSvc, "MucCalibMgr");
377 log << MSG::INFO << "Initialize NTuples" << endreq;
378
379 // Book ntuple
380 Gaudi::svcLocator() -> service("NTupleSvc", ntupleSvc);
381
382 StatusCode sc;
383
384 // event log
385 NTuplePtr nt1(ntupleSvc, "FILE450/EventLog");
386 if ( nt1 ) { m_eventLogTuple = nt1; }
387 else
388 {
389 m_eventLogTuple = ntupleSvc->book ("FILE450/EventLog", CLID_RowWiseTuple, "MucCalibConst N-Tuple");
390 if ( m_eventLogTuple )
391 {
392 sc = m_eventLogTuple->addItem ("event_id", m_ntEventId);
393 sc = m_eventLogTuple->addItem ("event_tag", m_ntEventTag);
394 sc = m_eventLogTuple->addItem ("start_time", m_ntEsTime);
395 sc = m_eventLogTuple->addItem ("digi_num", m_ntDigiNum);
396 sc = m_eventLogTuple->addItem ("track_num", m_ntTrackNum);
397 sc = m_eventLogTuple->addItem ("exphit_num", m_ntExpHitNum);
398 sc = m_eventLogTuple->addItem ("effhit_num", m_ntEffHitNum);
399 sc = m_eventLogTuple->addItem ("noshit_num", m_ntNosHitNum);
400 sc = m_eventLogTuple->addItem ("cluster_num", m_ntClusterNum);
401 sc = m_eventLogTuple->addItem ("event_time", m_ntEventTime);
402 }
403 else {
404 log << MSG::ERROR << "Cannot book N-tuple:" << long(m_eventLogTuple) << endmsg;
405 return StatusCode::FAILURE;
406 }
407 }
408
409 // track info
410
411 NTuplePtr nt2(ntupleSvc, "FILE450/MdcTrkInfo");
412 if ( nt2 ) { m_mdcTrkInfoTuple = nt2; }
413 else
414 {
415 m_mdcTrkInfoTuple = ntupleSvc->book ("FILE450/MdcTrkInfo", CLID_RowWiseTuple, "MucCalibConst N-Tuple");
416 if ( m_mdcTrkInfoTuple )
417 {
418 sc = m_mdcTrkInfoTuple->addItem ("charge", m_charge);
419 sc = m_mdcTrkInfoTuple->addItem ("mdcpx", m_mdcpx);
420 sc = m_mdcTrkInfoTuple->addItem ("mdcpy", m_mdcpy);
421 sc = m_mdcTrkInfoTuple->addItem ("mdcpz", m_mdcpz);
422 sc = m_mdcTrkInfoTuple->addItem ("mdcpt", m_mdcpt);
423 sc = m_mdcTrkInfoTuple->addItem ("mdcpp", m_mdcpp);
424 sc = m_mdcTrkInfoTuple->addItem ("mdcphi", m_mdcphi);
425 sc = m_mdcTrkInfoTuple->addItem ("mdctheta", m_mdctheta);
426 }
427 else {
428 log << MSG::ERROR << "Cannot book N-tuple:" << long(m_mdcTrkInfoTuple) << endmsg;
429 return StatusCode::FAILURE;
430 }
431 }
432
433
434 NTuplePtr nt3(ntupleSvc, "FILE450/TrackInfo");
435 if ( nt3 ) { m_trackInfoTuple = nt3; }
436 else
437 {
438 m_trackInfoTuple = ntupleSvc->book ("FILE450/TrackInfo", CLID_RowWiseTuple, "MucCalibConst N-Tuple");
439 if ( m_trackInfoTuple )
440 {
441 sc = m_trackInfoTuple->addItem ("track_event", m_ntTrackEvent);
442 sc = m_trackInfoTuple->addItem ("track_tag", m_ntTrackTag);
443 sc = m_trackInfoTuple->addItem ("track_hits", m_ntTrackHits);
444 sc = m_trackInfoTuple->addItem ("segment_fly", m_ntTrackSegFly);
445 sc = m_trackInfoTuple->addItem ("layer_fly_a", m_ntTrackLayFlyA);
446 sc = m_trackInfoTuple->addItem ("layer_fly_b", m_ntTrackLayFlyB);
447 sc = m_trackInfoTuple->addItem ("layer_fly_c", m_ntTrackLayFlyC);
448 sc = m_trackInfoTuple->addItem ("rec_mode", m_trkRecMode);
449 sc = m_trackInfoTuple->addItem ("chi2", m_chi2);
450 sc = m_trackInfoTuple->addItem ("px", m_px);
451 sc = m_trackInfoTuple->addItem ("py", m_py);
452 sc = m_trackInfoTuple->addItem ("pz", m_pz);
453 sc = m_trackInfoTuple->addItem ("pt", m_pt);
454 sc = m_trackInfoTuple->addItem ("pp", m_pp);
455 sc = m_trackInfoTuple->addItem ("r", m_r );
456 sc = m_trackInfoTuple->addItem ("costheta", m_cosTheta);
457 sc = m_trackInfoTuple->addItem ("theta", m_theta);
458 sc = m_trackInfoTuple->addItem ("phi", m_phi);
459 sc = m_trackInfoTuple->addItem ("depth", m_depth);
460 sc = m_trackInfoTuple->addItem ("br_last_lay", m_brLastLayer);
461 sc = m_trackInfoTuple->addItem ("ec_last_lay", m_ecLastLayer);
462 sc = m_trackInfoTuple->addItem ("total_hits", m_totalHits);
463 sc = m_trackInfoTuple->addItem ("fired_layers", m_totalLayers);
464 sc = m_trackInfoTuple->addItem ("maxhits_in_layer", m_maxHitsInLayer);
465 }
466 else {
467 log << MSG::ERROR << "Cannot book N-tuple:" << long(m_trackInfoTuple) << endmsg;
468 return StatusCode::FAILURE;
469 }
470 }
471
472 // track collinearity
473 NTuplePtr nt4(ntupleSvc, "FILE450/TrackDiff");
474 if ( nt4 ) { m_trackDiffTuple = nt4; }
475 else
476 {
477 m_trackDiffTuple = ntupleSvc->book ("FILE450/TrackDiff", CLID_RowWiseTuple, "MucCalibConst N-Tuple");
478 if ( m_trackDiffTuple ) {
479 sc = m_trackDiffTuple->addItem ("dimu_tag", m_ntDimuTag);
480 sc = m_trackDiffTuple->addItem ("pos_phi_diff", m_ntPosPhiDiff);
481 sc = m_trackDiffTuple->addItem ("pos_theta_diff", m_ntPosThetaDiff);
482 sc = m_trackDiffTuple->addItem ("mom_phi_diff", m_ntMomPhiDiff);
483 sc = m_trackDiffTuple->addItem ("mom_theta_diff", m_ntMomThetaDiff);
484 }
485 else {
486 log << MSG::ERROR << "Cannot book N-tuple:" << long(m_trackDiffTuple) << endmsg;
487 return StatusCode::FAILURE;
488 }
489 }
490
491 // cluster size
492 NTuplePtr nt5(ntupleSvc, "FILE450/ClusterSize");
493 if ( nt5 ) { m_clusterSizeTuple = nt5; }
494 else
495 {
496 m_clusterSizeTuple = ntupleSvc->book ("FILE450/ClusterSize", CLID_RowWiseTuple, "MucCalibConst N-Tuple");
497 if ( m_clusterSizeTuple ) {
498 sc = m_clusterSizeTuple->addItem ("cluster_size", m_ntClusterSize);
499 }
500 else {
501 log << MSG::ERROR << "Cannot book N-tuple:" << long(m_clusterSizeTuple) << endmsg;
502 return StatusCode::FAILURE;
503 }
504 }
505
506 // eff window
507 NTuplePtr nt6(ntupleSvc, "FILE450/EffWindow");
508 if ( nt6 ) { m_effWindowTuple = nt6; }
509 else
510 {
511 m_effWindowTuple = ntupleSvc->book ("FILE450/EffWindow", CLID_RowWiseTuple, "MucCalibConst N-Tuple");
512 if ( m_effWindowTuple ) {
513 sc = m_effWindowTuple->addItem ("hit_window", m_ntEffWindow);
514 }
515 else {
516 log << MSG::ERROR << "Cannot book N-tuple:" << long(m_effWindowTuple) << endmsg;
517 return StatusCode::FAILURE;
518 }
519 }
520 // res info
521 NTuplePtr nt7(ntupleSvc, "FILE450/ResInfo");
522 if ( nt7 ) { m_resInfoTuple = nt7; }
523 else
524 {
525 m_resInfoTuple = ntupleSvc->book ("FILE450/ResInfo", CLID_RowWiseTuple, "MucCalibConst N-Tuple");
526 if ( m_resInfoTuple ) {
527 sc = m_resInfoTuple->addItem ("line_res", m_lineRes);
528 sc = m_resInfoTuple->addItem ("quad_res", m_quadRes);
529 sc = m_resInfoTuple->addItem ("extr_res", m_extrRes);
530 sc = m_resInfoTuple->addItem ("res_part", m_resPart);
531 sc = m_resInfoTuple->addItem ("res_segment", m_resSegment);
532 sc = m_resInfoTuple->addItem ("res_layer", m_resLayer);
533 sc = m_resInfoTuple->addItem ("res_fired", m_resFired);
534 sc = m_resInfoTuple->addItem ("res_mode", m_resMode);
535 }
536 else {
537 log << MSG::ERROR << "Cannot book N-tuple:" << long(m_resInfoTuple) << endmsg;
538 return StatusCode::FAILURE;
539 }
540 }
541
542/*
543 NTuplePtr nt7(ntupleSvc, "FILE450/ResInfo");
544 if ( nt7 ) { m_resInfoTuple = nt7; }
545 else
546 {
547 m_resInfoTuple = ntupleSvc->book ("FILE450/ResInfo", CLID_ColumnWiseTuple, "MucCalibConst N-Tuple");
548 if ( m_resInfoTuple ) {
549 sc = m_resInfoTuple->addItem ("exp_num", m_nExpNum, 0, 30);
550 sc = m_resInfoTuple->addIndexedItem ("res", m_nExpNum, m_res);
551 sc = m_resInfoTuple->addIndexedItem ("res_part", m_nExpNum, m_resPart);
552 sc = m_resInfoTuple->addIndexedItem ("res_segment", m_nExpNum, m_resSegment);
553 sc = m_resInfoTuple->addIndexedItem ("res_layer", m_nExpNum, m_resLayer);
554 sc = m_resInfoTuple->addIndexedItem ("res_fired", m_nExpNum, m_resFired);
555 }
556 else {
557 log << MSG::ERROR << "Cannot book N-tuple:" << long(m_resInfoTuple) << endmsg;
558 return StatusCode::FAILURE;
559 }
560 }
561
562 for(int i=0; i<24; i++ )
563 {
564 m_res[i] = 211;
565 m_resPart[i] = -1;
566 m_resSegment[i] = -1;
567 m_resLayer[i] = -1;
568 m_resFired[i] = false;
569 }
570*/
571
572 return StatusCode::SUCCESS;
573}
574
575// Initialize histogram maps of fired strips and reconstructed hits
576StatusCode MucCalibMgr::InitHisto()
577{
578 MsgStream log(msgSvc, "MucCalibMgr");
579 log << MSG::INFO << "Initialize Histograms" << endreq;
580
581 // Init online monitor histo
582 InitOnlineHisto();
583
584 // Init eff map and so on
585 InitHistoLV2();
586 InitHistoLV1();
587 InitHistoLV0();
588
589 // Init 2D eff map
590 if( m_usePad != 0 ) Init2DEffHisto();
591
592 // Init cluster histo
593 InitClusterHisto();
594
595 // Init spacial resolution histo
596 InitResHisto();
597
598 return StatusCode::SUCCESS;
599}
600
601// Init histograme for online monitor
602StatusCode MucCalibMgr::InitOnlineHisto()
603{
604 char name[60];
605 char title[60];
606 int stripMax = 0;
607
608 // Init hit map
609 for( int i=0; i<B_LAY_NUM; i++ )
610 {
611 // According to layer
612 if(i%2 != 0 ) { stripMax = 96*7 + 112; }
613 else { stripMax = 48*8; }
614 sprintf( name, "HitMapBarrel_Lay_L%d", i );
615 sprintf( title, "Hit map in barrel layer %d", i );
616 m_hHitMapBarrel_Lay[i] = new TH1F(name,title, stripMax, 0, stripMax);
617
618 if( i<E_LAY_NUM )
619 {
620 stripMax = 64*4; // 256
621 sprintf( name, "HitMapEastEndcap_L%d", i );
622 sprintf( title, "Hit map in east-endcap layer %d", i );
623 m_hHitMapEndcap_Lay[0][i] = new TH1F(name,title, stripMax, 0, stripMax);
624
625 sprintf( name, "HitMapWestEndcap_L%d", i );
626 sprintf( title, "Hit map in west-endcap layer %d", i );
627 m_hHitMapEndcap_Lay[1][i] = new TH1F(name,title, stripMax, 0, stripMax);
628 }
629
630 }
631
632 for( int i=0; i<B_SEG_NUM; i++ )
633 {
634 // According to segment
635 sprintf( name, "HitMapBarrel_Seg_S%d", i );
636 sprintf( title, "Hit map in barrel segment %d", i );
637 if(i==2 ) { stripMax = 48*5 + 112*4; } // 688
638 else { stripMax = 48*5 + 96*4; } // 624
639 m_hHitMapBarrel_Seg[i] = new TH1F(name,title, stripMax, 0, stripMax);
640
641 if( i<E_SEG_NUM )
642 {
643 sprintf( name, "HitMapEastEndcap_S%d", i );
644 sprintf( title, "Hit map in east-endcap segment %d", i );
645 stripMax = 64*8; // 512
646 m_hHitMapEndcap_Seg[0][i] = new TH1F(name,title, stripMax, 0, stripMax);
647
648 sprintf( name, "HitMapWestEndcap_S%d", i );
649 sprintf( title, "Hit map in west-endcap segment %d", i );
650 m_hHitMapEndcap_Seg[1][i] = new TH1F(name,title, stripMax, 0, stripMax);
651 }
652 }
653
654 // Init histograms for online monitor
655 // 1D histograms
656 m_hHitVsEvent = new TH1F("HitVsEvent","Hit VS event",10000,0,10000);
657 m_hHitVsEvent->GetXaxis()->SetTitle("Event NO.");
658 m_hHitVsEvent->GetXaxis()->CenterTitle();
659 m_hHitVsEvent->GetYaxis()->SetTitle("Hits");
660 m_hHitVsEvent->GetYaxis()->CenterTitle();
661
662 m_hTrackDistance = new TH1F("TrackDistance","Track distance", 1000,-500,500);
663 m_hTrackDistance->GetXaxis()->SetTitle("Distance of fit pos and hit pos on 1st layer [cm]");
664 m_hTrackDistance->GetXaxis()->CenterTitle();
665
666 m_hTrackPosPhiDiff = new TH1F("TrackPosPhiDiff","#Delta#phi of tracks pos", 720,-2*PI,2*PI);
667 m_hTrackPosPhiDiff->GetXaxis()->SetTitle("#Delta#phi [rad]");
668 m_hTrackPosPhiDiff->GetXaxis()->CenterTitle();
669
670 m_hTrackPosThetaDiff = new TH1F("TrackPosThetaDiff","#Delta#theta of tracks pos", 720,-2*PI,2*PI);
671 m_hTrackPosThetaDiff->GetXaxis()->SetTitle("#Delta#theta [rad]");
672 m_hTrackPosThetaDiff->GetXaxis()->CenterTitle();
673
674 m_hTrackMomPhiDiff = new TH1F("TrackMomPhiDiff","#Delta#phi of tracks mom", 720,-2*PI,2*PI);
675 m_hTrackMomPhiDiff->GetXaxis()->SetTitle("#Delta#phi [rad]");
676 m_hTrackMomPhiDiff->GetXaxis()->CenterTitle();
677
678 m_hTrackMomThetaDiff = new TH1F("TrackMomThetaDiff","#Delta#theta of tracks mom", 720,-2*PI,2*PI);
679 m_hTrackMomThetaDiff->GetXaxis()->SetTitle("#Delta#theta [rad]");
680 m_hTrackMomThetaDiff->GetXaxis()->CenterTitle();
681
682 // 2D histograms
683 m_hDimuTracksPosDiff = new TH2F("DimuTracksPosDiff", "#Delta#phi VS #Delta#theta of dimu tracks pos", 720, -PI, PI, 720, -PI, PI);
684 m_hDimuTracksPosDiff->GetXaxis()->SetTitle("#Delta#theta");
685 m_hDimuTracksPosDiff->GetXaxis()->CenterTitle();
686 m_hDimuTracksPosDiff->GetYaxis()->SetTitle("#Delta#phi");
687 m_hDimuTracksPosDiff->GetYaxis()->CenterTitle();
688
689 m_hDimuTracksMomDiff = new TH2F("DimuTracksMomDiff", "#Delta#phi VS #Delta#theta of dimu tracks mom", 720, -PI, PI, 720, -PI, PI);
690 m_hDimuTracksMomDiff->GetXaxis()->SetTitle("#Delta#theta");
691 m_hDimuTracksMomDiff->GetXaxis()->CenterTitle();
692 m_hDimuTracksMomDiff->GetYaxis()->SetTitle("#Delta#phi");
693 m_hDimuTracksMomDiff->GetYaxis()->CenterTitle();
694
695 m_hPhiCosTheta = new TH2F("PhiVsCosTheta", "#phi VS cos(#theta)", 720, -1, 1, 720, -PI, PI);
696 m_hPhiCosTheta->GetXaxis()->SetTitle("cos(#theta)");
697 m_hPhiCosTheta->GetXaxis()->CenterTitle();
698 m_hPhiCosTheta->GetYaxis()->SetTitle("#phi");
699 m_hPhiCosTheta->GetYaxis()->CenterTitle();
700
701 return StatusCode::SUCCESS;
702}
703
704// Init histograme for calibration level 0
705StatusCode MucCalibMgr::InitHistoLV0()
706{
707 m_hBrLayerFire = new TH1F("BrLayerFire","Fires per layer in Barrel", LAYER_MAX+1, -0.5, LAYER_MAX+0.5);
708 m_hEcLayerFire = new TH1F("EcLayerFire","Fires per layer in Endcap", 2*LAYER_MAX-1, -LAYER_MAX+1, LAYER_MAX-0.5);
709
710 m_hLayerFire = new TH1F("LayerFire","Fires per layer", LAYER_MAX+1, -0.5, LAYER_MAX+0.5);
711 m_hLayerExpHit = new TH1F("LayerExpHit","Exp hits per layer", LAYER_MAX+1, -0.5, LAYER_MAX+0.5);
712 m_hLayerEffHit = new TH1F("LayerEffHit","Eff hits per layer", LAYER_MAX+1, -0.5, LAYER_MAX+0.5);
713 m_hLayerNosHit = new TH1F("LayerNosHit","Nos hits per layer", LAYER_MAX+1, -0.5, LAYER_MAX+0.5);
714 m_hLayerEff = new TH1F("LayerEff","Layer efficiency", LAYER_MAX+1, -0.5, LAYER_MAX+0.5);
715 m_hLayerNosRatio= new TH1F("LayerNosRatio","Layer noise hit ratio", LAYER_MAX+1, -0.5, LAYER_MAX+0.5);
716 m_hLayerArea = new TH1F("LayerArea","Layer area", LAYER_MAX+1, -0.5, LAYER_MAX+0.5);
717 m_hLayerNos = new TH1F("LayerNos","Layer noise", LAYER_MAX+1, -0.5, LAYER_MAX+0.5);
718 m_hLayerCnt = new TH1F("LayerCnt","Layer count", LAYER_MAX+1, -0.5, LAYER_MAX+0.5);
719
720 return StatusCode::SUCCESS;
721}
722
723// Init histograme for calibration level 1
724StatusCode MucCalibMgr::InitHistoLV1()
725{
726 m_hBoxFire = new TH1F("BoxFire","Fires per box", BOX_MAX+1, -0.5, BOX_MAX+0.5);
727 m_hBoxExpHit = new TH1F("BoxExpHit","Exp hits per box", BOX_MAX+1, -0.5, BOX_MAX+0.5);
728 m_hBoxEffHit = new TH1F("BoxEffHit","Eff hits per box", BOX_MAX+1, -0.5, BOX_MAX+0.5);
729 m_hBoxNosHit = new TH1F("BoxNosHit","Nos hits per box", BOX_MAX+1, -0.5, BOX_MAX+0.5);
730 m_hBoxEff = new TH1F("BoxEff","Box efficiency", BOX_MAX+1, -0.5, BOX_MAX+0.5);
731 m_hBoxNosRatio = new TH1F("BoxNosRatio","Box noise hit ratio", BOX_MAX+1, -0.5, BOX_MAX+0.5);
732 m_hBoxArea = new TH1F("BoxArea","Box area", BOX_MAX+1, -0.5, BOX_MAX+0.5);
733 m_hBoxNos = new TH1F("BoxNos","Box noise", BOX_MAX+1, -0.5, BOX_MAX+0.5);
734 m_hBoxCnt = new TH1F("BoxCnt","Box count", BOX_MAX+1, -0.5, BOX_MAX+0.5);
735
736 return StatusCode::SUCCESS;
737}
738
739// Init histograme for calibration level 2
740StatusCode MucCalibMgr::InitHistoLV2()
741{
742 char name[60];
743 char title[60];
744 int part, segment, layer, stripMax;
745 part = segment = layer = stripMax = 0;
746
747 for( int i=0; i<BOX_MAX; i++ )
748 {
749 m_ptrIdTr->SetBoxPos( i, &part, &segment, &layer );
750 stripMax = m_ptrIdTr->GetStripMax( part, segment, layer );
751
752 sprintf( name, "StripFireMap_P%d_S%d_L%d_Box%d", part, segment, layer, i );
753 sprintf( title, "Fires per strip in P%d_S%d_L%d Box%d",part, segment, layer, i );
754 m_hStripFireMap[i] = new TH1F(name,title, stripMax, 0, stripMax);
755
756 sprintf( name, "StripExpHitMap_P%d_S%d_L%d_Box%d", part, segment, layer, i );
757 sprintf( title, "Exp hits per strip in P%d_S%d_L%d Box%d",part, segment, layer, i );
758 m_hStripExpHitMap[i] = new TH1F(name,title, stripMax, 0, stripMax);
759
760 sprintf( name, "StripEffHitMap_P%d_S%d_L%d_Box%d", part, segment, layer, i );
761 sprintf( title, "Eff hits per strip in P%d_S%d_L%d Box%d",part, segment, layer, i );
762 m_hStripEffHitMap[i] = new TH1F(name,title, stripMax, 0, stripMax);
763
764 sprintf( name, "StripNosHitMap_P%d_S%d_L%d_Box%d", part, segment, layer, i );
765 sprintf( title, "Inc hits per strip in P%d_S%d_L%d Box%d",part, segment, layer, i );
766 m_hStripNosHitMap[i] = new TH1F(name,title, stripMax, 0, stripMax);
767
768 sprintf( name, "StripEffMap_P%d_S%d_L%d_Box%d", part, segment, layer, i );
769 sprintf( title, "Strip efficiency in P%d_S%d_L%d Box%d",part, segment, layer, i );
770 m_hStripEffMap[i] = new TH1F(name,title, stripMax, 0, stripMax);
771
772 sprintf( name, "StripNosRatioMap_P%d_S%d_L%d_Box%d", part, segment, layer, i );
773 sprintf( title, "Strip noise hit ratio in P%d_S%d_L%d Box%d",part, segment, layer, i );
774 m_hStripNosRatioMap[i] = new TH1F(name,title, stripMax, 0, stripMax);
775 }
776
777 m_hStripFire = new TH1F("StripFire", "Fires per strip", STRIP_MAX+1, -0.5, STRIP_MAX+0.5 );
778 m_hStripExpHit = new TH1F("StripExpHit", "Exp hit per strip", STRIP_MAX+1, -0.5, STRIP_MAX+0.5 );
779 m_hStripEffHit = new TH1F("StripEffHit", "Eff hit per strip", STRIP_MAX+1, -0.5, STRIP_MAX+0.5 );
780 m_hStripNosHit = new TH1F("StripNoshit", "Nos hit per strip", STRIP_MAX+1, -0.5, STRIP_MAX+0.5 );
781 m_hStripEff = new TH1F("StripEff", "Strip efficiency", STRIP_MAX+1, -0.5, STRIP_MAX+0.5 );
782 m_hStripNosRatio= new TH1F("StripNosRatio", "Strip noise hit ratio", STRIP_MAX+1, -0.5, STRIP_MAX+0.5 );
783 m_hStripArea = new TH1F("StripArea", "Strip area", STRIP_MAX+1, -0.5, STRIP_MAX+0.5 );
784 m_hStripNos = new TH1F("StripNos", "Strip noise", STRIP_MAX+1, -0.5, STRIP_MAX+0.5 );
785 m_hStripCnt = new TH1F("StripCnt", "Strip count", STRIP_MAX+1, -0.5, STRIP_MAX+0.5 );
786
787 return StatusCode::SUCCESS;
788}
789
790// Init 2D eff histograme
791StatusCode MucCalibMgr::Init2DEffHisto()
792{
793 char name[60];
794 int stripMax, boxID, stripID, xBin, yBin;
795 stripMax = boxID = stripID = xBin = yBin = 0;
796 double xStart, xEnd, yStart, yEnd, sWidth;
797 xStart = xEnd = yStart = yEnd = sWidth = 0.;
798
799 m_histArray = new TObjArray();
800
801 for( int i=0; i<PART_MAX; i++ )
802 for( int j=0; j<((i==BRID)?B_SEG_NUM:E_SEG_NUM); j++ )
803 for( int k=0; k<((i==BRID)?B_LAY_NUM:E_LAY_NUM); k++ )
804 {
805 boxID = m_ptrIdTr->GetBoxId(i,j,k);
806
807 if( i==BRID )
808 {
809 xStart = -2000, xEnd = 2000;
810 sWidth = B_STR_DST[k];
811 xBin = int((xEnd - xStart)/sWidth);
812 yStart = -B_BOX_WT[k]/2-100, yEnd = B_BOX_WT[k]/2+100;
813 yBin = int((B_BOX_WT[k]+200)/sWidth);
814 }
815 else
816 {
817 xStart = yStart = -1250;
818 xEnd = yEnd = 1250;
819 sWidth = E_STR_DST;
820 xBin = yBin = int((xEnd - xStart)/sWidth);
821 }
822
823 sprintf(name, "ExpMap2D_P%d_S%d_L%d_Box%d", i, j, k, boxID);
824 m_h2DExpMap[i][j][k] = new TH2F(name, name, xBin, xStart, xEnd, yBin, yStart, yEnd);
825 sprintf(name, "HitMap2D_P%d_S%d_L%d_Box%d", i, j, k, boxID);
826 m_h2DHitMap[i][j][k] = new TH2F(name, name, xBin, xStart, xEnd, yBin, yStart, yEnd);
827 sprintf(name, "EffMap2D_P%d_S%d_L%d_Box%d", i, j, k, boxID);
828 m_h2DEffMap[i][j][k] = new TH2F(name, name, xBin, xStart, xEnd, yBin, yStart, yEnd);
829
830 //m_histArray->Add(m_h2DExpMap[i][j][k]);
831 //m_histArray->Add(m_h2DHitMap[i][j][k]);
832 m_histArray->Add(m_h2DEffMap[i][j][k]);
833 }
834
835 return StatusCode::SUCCESS;
836}
837
838// Init histograme for cluster
839StatusCode MucCalibMgr::InitClusterHisto()
840{
841 char name[60];
842 char title[60];
843 int part, segment, layer, stripMax;
844 part = segment = layer = stripMax = 0;
845
846 for( int i=0; i<LAYER_MAX; i++ )
847 {
848 sprintf( name, "LayerCluster_L%d", i );
849 sprintf( title, "Clusters in L%d",i );
850 m_hLayerCluster[i] = new TH1F(name,title, CLUSTER_RANGE, 0, CLUSTER_RANGE );
851 }
852
853 for( int i=0; i<BOX_MAX; i++ )
854 {
855 m_ptrIdTr->SetBoxPos( i, &part, &segment, &layer );
856 stripMax = m_ptrIdTr->GetStripMax( part, segment, layer );
857 sprintf( name, "BoxCluster_P%d_S%d_L%d_Box%d", part, segment, layer, i );
858 sprintf( title, "Clusters in P%d_S%d_L%d Box%d", part, segment, layer, i );
859 m_hBoxCluster[i] = new TH1F(name,title, CLUSTER_RANGE, 0, CLUSTER_RANGE );
860 }
861
862 m_hLayerClusterCmp = new TH1F("LayerCluster", "LayerCluster", LAYER_MAX+1, -0.5, LAYER_MAX+0.5 );
863 m_hBoxClusterCmp = new TH1F("BoxCluster", "BoxCluster", BOX_MAX+1, -0.5, BOX_MAX+0.5 );
864
865 return StatusCode::SUCCESS;
866}
867
868// Init histograms for spacial resolution
869StatusCode MucCalibMgr::InitResHisto()
870{
871 char name[60];
872 char title[60];
873
874 for( int i=0; i<B_LAY_NUM; i++ )
875 {
876 sprintf( name, "BarrelRes_L%d", i );
877 sprintf( title, "Barrel spacial resolution in L%d",i );
878 m_hBarrelResDist[i] = new TH1F(name,title, 200, -100, 100 );
879 }
880
881 for( int i=0; i<E_LAY_NUM; i++ )
882 {
883 sprintf( name, "EndcapRes_L%d", i );
884 sprintf( title, "Endcap spacial resolution in L%d",i );
885 m_hEndcapResDist[i] = new TH1F(name,title, 200, -100, 100 );
886 }
887
888 m_hBarrelResComp[0] = new TH1F("BarrelResMean", "BarrelResMean", B_LAY_NUM+1, -0.5, B_LAY_NUM+0.5 );
889 m_hBarrelResComp[1] = new TH1F("BarrelResSigma", "BarrelResSigma", B_LAY_NUM+1, -0.5, B_LAY_NUM+0.5 );
890 m_hEndcapResComp[0] = new TH1F("EndcapResMean", "EndcapResMean", E_LAY_NUM+1, -0.5, E_LAY_NUM+0.5 );
891 m_hEndcapResComp[1] = new TH1F("EndcapResSigma", "EndcapResSigma", E_LAY_NUM+1, -0.5, E_LAY_NUM+0.5 );
892
893 return StatusCode::SUCCESS;
894}
895
896// Init Tree
897StatusCode MucCalibMgr::InitConstTree()
898{
899 MsgStream log(msgSvc, "MucCalibMgr");
900 log << MSG::INFO << "Initialize Trees" << endreq;
901
902 // job log tree
903 m_tJobLog = new TTree("JobLog","Job information");
904 m_tJobLog->Branch("version", &m_vs, "m_vs/D");
905 m_tJobLog->Branch("high_voltage", &m_hv, "m_hv/D");
906 m_tJobLog->Branch("threshold", &m_th, "m_th/D");
907 m_tJobLog->Branch("event_rate", &m_er, "m_er/D");
908 m_tJobLog->Branch("input_tag", &m_tag, "m_tag/D");
909 m_tJobLog->Branch("rec_mode", &m_recMode, "m_recMode/I");
910 m_tJobLog->Branch("use_pad", &m_usePad, "m_usePad/I");
911 m_tJobLog->Branch("eff_window", &m_effWindow, "m_effWindow/I");
912 m_tJobLog->Branch("cluster_mode", &m_clusterMode, "m_clusterMode/I");
913 m_tJobLog->Branch("check_event", &m_checkEvent, "m_checkEvent/I");
914 m_tJobLog->Branch("dimu_select", &m_dimuSelect, "m_dimuSelect/I");
915 m_tJobLog->Branch("dimu_only", &m_dimuOnly, "m_dimuOnly/I");
916
917 m_tJobLog->Branch("run_start", &m_fStartRun, "m_fStartRun/I");
918 m_tJobLog->Branch("run_end", &m_fEndRun, "m_fEndRun/I");
919 m_tJobLog->Branch("daq_time", &m_fTotalDAQTime, "m_fTotalDAQTime/D");
920 m_tJobLog->Branch("job_time", &m_fTotalJobTime, "m_fTotalJobTime/D");
921 m_tJobLog->Branch("calib_layer", &m_fCalibLayerNum, "m_fCalibLayerNum/D");
922 m_tJobLog->Branch("calib_box", &m_fCalibBoxNum, "m_fCalibBoxNum/D");
923 m_tJobLog->Branch("calib_strip", &m_fCalibStripNum, "m_fCalibStripNum/D");
924 m_tJobLog->Branch("total_event", &m_fTotalEvent, "m_fTotalEvent/D");
925 m_tJobLog->Branch("total_digi", &m_fTotalDigi, "m_fTotalDigi/D");
926 m_tJobLog->Branch("total_exphit", &m_fTotalExpHit, "m_fTotalExpHit/D");
927 m_tJobLog->Branch("total_effhit", &m_fTotalEffHit, "m_fTotalEffHit/D");
928 m_tJobLog->Branch("total_noshit", &m_fTotalNosHit, "m_fTotalNosHit/D");
929 m_tJobLog->Branch("total_cluster", &m_fTotalClstNum, "m_fTotalClstNum/D");
930 m_tJobLog->Branch("total_strip_area",&m_fTotalStripArea,"m_fTotalStripArea/D");
931 m_tJobLog->Branch("layer_coverage", &m_fLayerCoverage, "m_fLayerCoverage/D");
932 m_tJobLog->Branch("box_coverage", &m_fBoxCoverage, "m_fBoxCoverage/D");
933 m_tJobLog->Branch("strip_coverage", &m_fStripCoverage, "m_fStripCoverage/D");
934
935 // statistic log tree
936 m_tStatLog = new TTree("StatLog","Statistic results");
937
938 // layer constants tree, level 0
939 m_tLayConst = new TTree("LayConst","layer constants");
940 m_tLayConst->Branch("layer_id", &m_fLayerId, "m_fLayerId/D");
941 m_tLayConst->Branch("layer_eff", &m_fLayerEff, "m_fLayerEff/D");
942 m_tLayConst->Branch("layer_efferr", &m_fLayerEffErr, "m_fLayerEffErr/D");
943 m_tLayConst->Branch("layer_nosratio", &m_fLayerNosRatio,"m_fLayerNosRatio/D");
944 m_tLayConst->Branch("layer_digi", &m_fLayerDigi, "m_fLayerDigi/D");
945 m_tLayConst->Branch("layer_noise", &m_fLayerNos, "m_fLayerNos/D");
946 m_tLayConst->Branch("layer_cnt", &m_fLayerCnt, "m_fLayerCnt/D");
947 m_tLayConst->Branch("layer_exphit", &m_fLayerExpHit, "m_fLayerExpHit/D");
948 m_tLayConst->Branch("layer_effHit", &m_fLayerEffHit, "m_fLayerEffHit/D");
949 m_tLayConst->Branch("layer_nosHit", &m_fLayerNosHit, "m_fLayerNosHit/D");
950 m_tLayConst->Branch("layer_cluster", &m_fLayerCluster, "m_fLayerCluster/D");
951 m_tLayConst->Branch("layer_trknum", &m_fLayerTrkNum, "m_fLayerTrkNum/D");
952
953 // box constants tree, level 1
954 m_tBoxConst = new TTree("BoxConst","box constants");
955 m_tBoxConst->Branch("box_id", &m_fBoxId, "m_fBoxId/D");
956 m_tBoxConst->Branch("box_part", &m_fBoxPart, "m_fBoxPart/D");
957 m_tBoxConst->Branch("box_segment", &m_fBoxSegment, "m_fBoxSegment/D");
958 m_tBoxConst->Branch("box_layer", &m_fBoxLayer, "m_fBoxLayer/D");
959 m_tBoxConst->Branch("box_eff", &m_fBoxEff, "m_fBoxEff/D");
960 m_tBoxConst->Branch("box_efferr", &m_fBoxEffErr, "m_fBoxEffErr/D");
961 m_tBoxConst->Branch("box_nosratio", &m_fBoxNosRatio,"m_fBoxNosRatio/D");
962 m_tBoxConst->Branch("box_digi", &m_fBoxDigi, "m_fBoxDigi/D");
963 m_tBoxConst->Branch("box_noise", &m_fBoxNos, "m_fBoxNos/D");
964 m_tBoxConst->Branch("box_cnt", &m_fBoxCnt, "m_fBoxCnt/D");
965 m_tBoxConst->Branch("box_exphit", &m_fBoxExpHit, "m_fBoxExpHit/D");
966 m_tBoxConst->Branch("box_effhit", &m_fBoxEffHit, "m_fBoxEffHit/D");
967 m_tBoxConst->Branch("box_noshit", &m_fBoxNosHit, "m_fBoxNosHit/D");
968 m_tBoxConst->Branch("box_cluster", &m_fBoxCluster, "m_fBoxCluster/D");
969 m_tBoxConst->Branch("box_trknum", &m_fBoxTrkNum, "m_fBoxTrkNum/D");
970
971 // strip constants tree, level 2
972 m_tStrConst = new TTree("StrConst","strip constants");
973 m_tStrConst->Branch("strip_id", &m_fStripId, "m_fStripId/D");
974 m_tStrConst->Branch("strip_part", &m_fStripPart, "m_fStripPart/D");
975 m_tStrConst->Branch("strip_segment", &m_fStripSegment, "m_fStripSegment/D");
976 m_tStrConst->Branch("strip_layer", &m_fStripLayer, "m_fStripLayer/D");
977 m_tStrConst->Branch("strip_eff", &m_fStripEff, "m_fStripEff/D");
978 m_tStrConst->Branch("strip_efferr", &m_fStripEffErr, "m_fStripEffErr/D");
979 m_tStrConst->Branch("strip_nosratio", &m_fStripNosRatio, "m_fStripNosRatio/D");
980 m_tStrConst->Branch("strip_digi", &m_fStripDigi, "m_fStripDigi/D");
981 m_tStrConst->Branch("strip_noise", &m_fStripNos, "m_fStripNos/D");
982 m_tStrConst->Branch("strip_cnt", &m_fStripCnt, "m_fStripCnt/D");
983 m_tStrConst->Branch("strip_effhit", &m_fStripEffHit, "m_fStripEffHit/D");
984 m_tStrConst->Branch("strip_exphit", &m_fStripExpHit, "m_fStripExpHit/D");
985 m_tStrConst->Branch("strip_noshit", &m_fStripNosHit, "m_fStripNosHit/D");
986 m_tStrConst->Branch("strip_trknum", &m_fStripTrkNum, "m_fStripTrkNum/D");
987
988 return StatusCode::SUCCESS;
989}
990
991// Init area of units
992StatusCode MucCalibMgr::InitArea()
993{
994 int stripMax, boxID, stripID;
995 stripMax = boxID = stripID = 0;
996 double totalArea = 0;
997
998 for( int i=0; i<PART_MAX; i++ )
999 for( int j=0; j<((i==BRID)?B_SEG_NUM:E_SEG_NUM); j++ )
1000 for( int k=0; k<((i==BRID)?B_LAY_NUM:E_LAY_NUM); k++ )
1001 {
1002 MucBoxCal* aBox = new MucBoxCal( i, j, k );
1003 boxID = m_ptrIdTr->GetBoxId(i, j, k);
1004 m_boxResults[3][boxID] = aBox->GetArea();
1005 m_layerResults[3][k] += aBox->GetArea();
1006 delete aBox;
1007
1008 stripMax = m_ptrIdTr->GetStripMax( i, j, k );
1009 for( int m=0; m<stripMax; m++ )
1010 {
1011 MucStripCal* aStrip = new MucStripCal( i, j, k, m );
1012 stripID =m_ptrIdTr->GetStripId(i, j, k, m );
1013 m_stripResults[3][stripID] = aStrip->GetArea();
1014 totalArea += m_stripResults[3][stripID];
1015 delete aStrip;
1016 }
1017 }
1018
1019 for( int i=0; i<LAYER_MAX; i++ ) m_hLayerArea->Fill(i, m_layerResults[3][i]);
1020 for( int i=0; i<BOX_MAX; i++ ) m_hBoxArea->Fill(i, m_boxResults[3][i]);
1021 for( int i=0; i<STRIP_MAX; i++ ) m_hStripArea->Fill(i, m_stripResults[3][i]);
1022
1023 m_fTotalStripArea = totalArea;
1024
1025 return StatusCode::SUCCESS;
1026}
1027
1028//-------------------------------------------------------------------------------------------------
1029//-------------------------------- Member functions for executing --------------------------------
1030//-------------------------------------------------------------------------------------------------
1031
1032// Dimu select
1033StatusCode MucCalibMgr::DimuSelect()
1034{
1035 MsgStream log(msgSvc, "MucCalibMgr");
1036 log << MSG::INFO << "Dimu select" << endreq;
1037 Gaudi::svcLocator()->service("EventDataSvc", eventSvc);
1038
1039 if( m_tag > 0) { m_eventTag = (int)m_tag; return ( StatusCode::SUCCESS ); }
1040
1041 m_eventTag = 0;
1042 SmartDataPtr<Event::EventHeader> eventHeader(eventSvc,"/Event/EventHeader");
1043 if(!eventHeader) {
1044 log << MSG::FATAL << "Could not find event header" << endreq;
1045 return( StatusCode::FAILURE );
1046 }
1047
1048 // Select by MDC Info
1049 SmartDataPtr<RecMdcTrackCol> mdcTrackCol(eventSvc,"/Event/Recon/RecMdcTrackCol");
1050 if(!mdcTrackCol) {
1051 log << MSG::FATAL << "Could not find mdc tracks" << endreq;
1052 return( StatusCode::FAILURE);
1053 }
1054 log << MSG::INFO << mdcTrackCol->size() << endreq;
1055 if( mdcTrackCol->size() != 2 ) return ( StatusCode::FAILURE );
1056
1057 log << MSG::INFO << "r1:\t" << (*mdcTrackCol)[0]->r() << "\tz1:" << (*mdcTrackCol)[0]->z() << endreq;
1058 log << MSG::INFO << "r2:\t" << (*mdcTrackCol)[1]->r() << "\tz2:" << (*mdcTrackCol)[1]->z() << endreq;
1059 log << MSG::INFO << "p1:\t" << (*mdcTrackCol)[0]->p() << "\tp2:" << (*mdcTrackCol)[1]->p() << endreq;
1060
1061 bool mdcFlag1 = (*mdcTrackCol)[0]->charge() + (*mdcTrackCol)[1]->charge() == 0;
1062 bool mdcFlag2 = fabs((*mdcTrackCol)[0]->r())<=1 && fabs((*mdcTrackCol)[0]->z())<3;
1063 bool mdcFlag3 = fabs((*mdcTrackCol)[1]->r())<=1 && fabs((*mdcTrackCol)[1]->z())<3;
1064 bool mdcFlag4 = (*mdcTrackCol)[0]->p()<2 && (*mdcTrackCol)[1]->p()<2;
1065 bool mdcFlag5 = fabs( (*mdcTrackCol)[0]->p() - (*mdcTrackCol)[1]->p() )/( (*mdcTrackCol)[0]->p() + (*mdcTrackCol)[1]->p() ) < 0.5;
1066 bool mdcPass = false;
1067 if( mdcFlag1 && mdcFlag2 && mdcFlag3 && mdcFlag4 && mdcFlag5 ) mdcPass = true;
1068
1069 // Select by TOF Info
1070 SmartDataPtr<RecTofTrackCol> tofTrackCol(eventSvc,"/Event/Recon/RecTofTrackCol");
1071 if(!tofTrackCol) {
1072 log << MSG::FATAL << "Could not find RecTofTrackCol!" << endreq;
1073 return( StatusCode::FAILURE);
1074 }
1075 log << MSG::INFO << "TOF tracks:\t" << tofTrackCol->size() << endreq;
1076
1077 double t1, t2;
1078 t1 = 0., t2 = 1000;
1079 // if(tofTrackCol->size() < 8 && tofTrackCol->size() > 20)
1080 if(tofTrackCol->size() > 7 && tofTrackCol->size() < 21)
1081 {
1082 int goodtof = 0;
1083 for(unsigned int itof = 0; itof < tofTrackCol->size(); itof++)
1084 {
1085 TofHitStatus *status = new TofHitStatus;
1086 status->setStatus((*tofTrackCol)[itof]->status());
1087
1088 if( !(status->is_cluster()) ) { delete status; continue; }
1089 if(goodtof==0) t1 = (*tofTrackCol)[itof]->tof();
1090 if(goodtof==1) t2 = (*tofTrackCol)[itof]->tof();
1091
1092 goodtof++;
1093 delete status;
1094 }
1095 }
1096 bool tofPass = false;
1097 if( fabs( t1-t2 ) < 4) tofPass = true; // ns
1098
1099 // Select by EMC Info
1100 SmartDataPtr<RecEmcShowerCol> emcShowerCol(eventSvc,"/Event/Recon/RecEmcShowerCol");
1101 if (!emcShowerCol) {
1102 log << MSG::FATAL << "Could not find RecEmcShowerCol!" << endreq;
1103 return( StatusCode::FAILURE);
1104 }
1105 log << MSG::INFO << "EMC showers:\t" << emcShowerCol->size() << endreq;
1106
1107 if( emcShowerCol->size() < 2 ) return ( StatusCode::SUCCESS );
1108
1109 double e1, e2, theta1, theta2, phi1, phi2;
1110 int part;
1111
1112 e1 = (*emcShowerCol)[0]->energy(); e2 = (*emcShowerCol)[1]->energy();
1113 theta1 = (*emcShowerCol)[0]->theta(); theta2 = (*emcShowerCol)[1]->theta();
1114 phi1 = (*emcShowerCol)[0]->phi(); phi2 = (*emcShowerCol)[1]->phi();
1115 part = (*emcShowerCol)[0]->module();
1116
1117 log << MSG::INFO << "e1:\t" << e1 << "\te2:\t" << e2 << endreq;
1118 log << MSG::INFO << "theta1:\t" << theta1 << "\ttheta2:\t" << theta2 << endreq;
1119 log << MSG::INFO << "phi1:\t" << phi1 << "\tphi2:\t" << phi2 << endreq;
1120
1121 bool emcFlag1 = e1 < 1.0 && e1 > 0.1;
1122 bool emcFlag2 = e2 < 1.0 && e2 > 0.1;
1123 bool emcFlag3 = fabs(theta1 + theta2 - PI) < 0.05;
1124 bool emcFlag4 = fabs(phi1 - phi2) - PI < 0.4;
1125
1126 bool emcPass = false;
1127 if( emcFlag1 && emcFlag2 && emcFlag3 && emcFlag4 ) emcPass = true;
1128
1129 // Select by MUC Info
1130 SmartDataPtr<MucDigiCol> mucDigiCol(eventSvc,"/Event/Digi/MucDigiCol");
1131 if(!mucDigiCol) {
1132 log << MSG::FATAL << "Could not find MUC digi" << endreq;
1133 return( StatusCode::FAILURE);
1134 }
1135 SmartDataPtr<RecMucTrackCol> mucTrackCol(eventSvc,"/Event/Recon/RecMucTrackCol");
1136 if (!mucTrackCol) {
1137 log << MSG::FATAL << "Could not find RecMucTrackCol" << endreq;
1138 return( StatusCode::FAILURE);
1139 }
1140 log << MSG::INFO << "digi:\t" << mucDigiCol->size() << "tracks:\t" << mucTrackCol->size() << endreq;
1141
1142 bool mucFlag1 = mucDigiCol->size()>6;
1143 bool mucFlag2 = mucTrackCol->size()>1;
1144 bool mucPass = false;
1145 if( mucFlag1 && mucFlag2 ) mucPass = true;
1146
1147 if( mdcPass && tofPass && emcPass && mucPass ) m_eventTag = 1;
1148 else m_eventTag = (int)m_tag;
1149
1150 return( StatusCode::SUCCESS );
1151}
1152
1153// Read event info from TDS
1154StatusCode MucCalibMgr::ReadEvent()
1155{
1156 MsgStream log(msgSvc, "MucCalibMgr");
1157 log << MSG::INFO << "Read event" << endreq;
1158
1159 Gaudi::svcLocator()->service("EventDataSvc", eventSvc);
1160 m_evtBegin = clock();
1161
1162 // Check event head
1163 SmartDataPtr<Event::EventHeader> eventHeader(eventSvc,"/Event/EventHeader");
1164 if(!eventHeader) {
1165 log << MSG::FATAL << "Could not find event header" << endreq;
1166 return( StatusCode::FAILURE );
1167 }
1168
1169 m_currentRun = eventHeader->runNumber();
1170 m_currentEvent = eventHeader->eventNumber();
1171 if( m_fStartRun == 0 ) m_fStartRun = m_currentRun;
1172 m_fEndRun = m_currentRun;
1173 m_fTotalEvent++;
1174
1175 log << MSG::INFO << "Run [ " << m_currentRun << " ]\tEvent [ " << m_currentEvent << " ]" << endreq;
1176 if( ((long)m_fTotalEvent)%2000 == 0 ) cout << m_fTotalEvent << "\tdone!" << endl;
1177
1178 // Select dimu
1179 if( m_dimuSelect ) {
1180 MucCalibMgr::DimuSelect();
1181 log << MSG::INFO << "Event tag:\t" << m_eventTag << endreq;
1182 if( m_dimuOnly && m_eventTag != 1 ) return( StatusCode::FAILURE );
1183 }
1184
1185 //---> Retrieve MUC digi
1186 log << MSG::INFO << "Retrieve digis" << endreq;
1187
1188 SmartDataPtr<MucDigiCol> mucDigiCol(eventSvc,"/Event/Digi/MucDigiCol");
1189 if(!mucDigiCol) {
1190 log << MSG::FATAL << "Could not find MUC digi" << endreq;
1191 return( StatusCode::FAILURE);
1192 }
1193
1194 int part, segment, layer, strip, pad;
1195 part = segment = layer = strip = pad = 0;
1196 double padX, padY, padZ;
1197 padX = padY = padZ = 0.;
1198 double resMax = 0.;
1199
1200 Identifier mucId;
1201 MucDigiCol::iterator digiIter = mucDigiCol->begin();
1202 int eventDigi = 0;
1203 for ( int digiId =0; digiIter != mucDigiCol->end(); digiIter++, digiId++ )
1204 {
1205 mucId = (*digiIter)->identify();
1206 part = MucID::barrel_ec(mucId);
1207 segment = MucID::segment(mucId);
1208 layer = MucID::layer(mucId);
1209 strip = MucID::channel(mucId);
1210
1211 log << MSG::DEBUG << "[" << part << "\t" << segment << "\t" << layer << "\t" << strip << "]\t" ;
1212 if( (digiId+1)%8 == 0 ) log << MSG::DEBUG << endreq;
1213
1214 eventDigi ++;
1215
1216 if(abs(part)>=PART_MAX || abs(segment)>=SEGMENT_MAX || abs(layer)>=LAYER_MAX || abs(strip)>=STRIP_INBOX_MAX) {
1217 log << MSG::ERROR << endreq << "Digi IDs slop over!" << endreq;
1218 continue;
1219 }
1220
1221 // Add digi
1222 MucMark *aMark = new MucMark( part, segment, layer, strip );
1223 m_digiCol.push_back( aMark );
1224 m_segDigiCol[part][segment].push_back( aMark );
1225 }
1226 log << MSG::DEBUG << endreq;
1227 log << MSG::INFO << "Total digits of this event: " << eventDigi << endreq;
1228 if( eventDigi > 200) log << MSG::ERROR << "Event: " << m_currentEvent << "\tdigits sharply rise:\t" << eventDigi << endreq;
1229
1230 /*
1231 if( (long)m_fTotalEvent%10000 == 0 ) {
1232 log << MSG::INFO << "Delete HitVsEvent." << endreq;
1233 if(m_hHitVsEvent != NULL) delete m_hHitVsEvent;
1234 m_hHitVsEvent = new TH1F("HitVsEvent","Hit VS event",10000,m_fTotalEvent,m_fTotalEvent+10000);
1235 m_hHitVsEvent->GetXaxis()->SetTitle("Event NO.");
1236 log << MSG::INFO << "Recreate HitVsEvent." << endreq;
1237 }
1238 log << MSG::INFO << "Fill HitVsEvent\t" << m_hHitVsEvent << "\t" << m_fTotalEvent << "\t" << m_currentEvent << endreq;
1239 //m_hHitVsEvent->Fill(m_currentEvent+m_currentEvent%10000, eventDigi);
1240 m_hHitVsEvent->Fill(m_fTotalEvent, eventDigi);
1241 log << MSG::INFO << "Fill HitVsEvent done." << endreq;
1242 */
1243
1244 // Search cluster in digis
1245 // Method detail in MucMark class
1246 int clusterNum, bigClusterNum, clusterSize;
1247 clusterNum = bigClusterNum = clusterSize = 0;
1248 if( m_clusterMode ) {
1249 log << MSG::INFO << "Searching clusters" << endreq;
1250 m_clusterCol = (*m_ptrMucMark).CreateClusterCol(m_clusterMode, m_digiCol );
1251 }
1252
1253 for( unsigned int i=0; i<m_clusterCol.size(); i++ )
1254 {
1255 clusterSize = m_clusterCol[i].size();
1256 // real cluster, size >= 2
1257 if( clusterSize > CLUSTER_ALARM )
1258 {
1259 log << MSG::WARNING << "Big cluster:" << endreq;
1260 part = (*m_clusterCol[i][0]).Part();
1261 segment = (*m_clusterCol[i][0]).Segment();
1262 layer = (*m_clusterCol[i][0]).Layer();
1263
1264 if( m_clusterSave ) (*m_fdata) << "Event:\t" << m_currentEvent << "\tbig cluster " << bigClusterNum << endl;
1265
1266 for( int j=0; j<clusterSize; j++ )
1267 {
1268 strip = (*m_clusterCol[i][j]).Strip();
1269 log << MSG::WARNING << "[" << part << "\t" << segment << "\t" << layer << "\t" << strip << "]\t";
1270 if( (j+1)%8 == 0 ) log << MSG::WARNING << endreq;
1271 if( m_clusterSave ) (*m_fdata) << part << "\t" << segment << "\t" << layer << "\t" << strip << endl;
1272 }
1273 log << MSG::WARNING << endreq;
1274 bigClusterNum ++;
1275 }
1276 else if( clusterSize > 1 )
1277 {
1278 log << MSG::DEBUG << "cluster: " << clusterNum << endreq;
1279 clusterNum ++, m_fTotalClstNum ++;
1280 part = (*m_clusterCol[i][0]).Part();
1281 segment = (*m_clusterCol[i][0]).Segment();
1282 layer = (*m_clusterCol[i][0]).Layer();
1283 for( int j=0; j<clusterSize; j++ )
1284 {
1285 strip = (*m_clusterCol[i][j]).Strip();
1286 log << MSG::DEBUG << "[" << part << "\t" << segment << "\t" << layer << "\t" << strip << "]\t";
1287 if( (j+1)%8 == 0 ) log << MSG::DEBUG << endreq;
1288 }
1289 log << MSG::DEBUG << endreq;
1290 }
1291 } // End m_clusterCol.size()
1292
1293 if( m_clusterMode) log << MSG::INFO << "Total clusters in this event: " << clusterNum << endreq;
1294 else log << MSG::INFO << "Clusters not built" << endreq;
1295 //<--- End retrieve digis
1296
1297 //---> Retrieve rec tracks
1298 log << MSG::INFO << "Retrieve tracks" << endreq;
1299 // MDC tracks
1300 SmartDataPtr<RecMdcTrackCol> mdcTrackCol(eventSvc,"/Event/Recon/RecMdcTrackCol");
1301 if(!mdcTrackCol) {
1302 log << MSG::FATAL << "Could not find mdc tracks" << endreq;
1303 return( StatusCode::FAILURE);
1304 }
1305
1306 RecMdcTrackCol::iterator mdctrkIter = mdcTrackCol->begin();
1307 for (; mdctrkIter != mdcTrackCol->end(); mdctrkIter++)
1308 {
1309 m_charge = (*mdctrkIter)->charge();
1310 m_mdcpx = (*mdctrkIter)->px();
1311 m_mdcpy = (*mdctrkIter)->py();
1312 m_mdcpz = (*mdctrkIter)->pz();
1313 m_mdcpt = (*mdctrkIter)->pxy();
1314 m_mdcpp = (*mdctrkIter)->p();
1315 m_mdcphi = (*mdctrkIter)->phi();
1316 m_mdctheta = (*mdctrkIter)->theta();
1317 m_mdcTrkInfoTuple->write();
1318 }
1319
1320 // MUC tracks
1321 SmartDataPtr<RecMucTrackCol> mucTrackCol(eventSvc,"/Event/Recon/RecMucTrackCol");
1322 if (!mucTrackCol) {
1323 log << MSG::FATAL << "Could not find RecMucTrackCol" << endreq;
1324 return( StatusCode::FAILURE);
1325 }
1326
1327 RecMucTrackCol *aRecMucTrackCol = mucTrackCol;
1328 if (aRecMucTrackCol->size() < 1) {
1329 log << MSG::INFO << "No MUC tracks in this event" << endreq;
1330 return StatusCode::SUCCESS;
1331 }
1332 log << MSG::INFO << "Total tracks of this event: " << aRecMucTrackCol->size() << endreq;
1333
1334 // Get RecEsTimeCol
1335 //int esTimeflag;
1336 //if( m_recMode == 0 ) // only for ExtTrk
1337 if( 0 ) // only for ExtTrk
1338 {
1339 SmartDataPtr<RecEsTimeCol> aRecEsTimeCol(eventSvc,"/Event/Recon/RecEsTimeCol");
1340 if( ! aRecEsTimeCol ){
1341 log << MSG::ERROR << "Could not find RecEsTimeCol" << endreq;
1342 return StatusCode::FAILURE;
1343 }else{
1344 RecEsTimeCol::iterator iter_evt = aRecEsTimeCol->begin();
1345 //tes = (*iter_evt)->getTest();
1346 //esTimeflag = (*iter_evt)->getStat();
1347 m_ntEsTime = (*iter_evt)->getStat();
1348 if( (*iter_evt)->getStat() != 211 ) {
1349 log << MSG::WARNING << "Event time not by TOF, skip!" << endreq;
1350 return StatusCode::SUCCESS;
1351 }
1352 }
1353 }
1354
1355 // Phi diff of two tracks( dimuon event )
1356 double phi1, phi2, phiDiff, theta1, theta2, thetaDiff;
1357 phiDiff = thetaDiff = 0.;
1358 if( aRecMucTrackCol->size()==2 && (*aRecMucTrackCol)[0]->GetTotalHits() > 4 && (*aRecMucTrackCol)[1]->GetTotalHits() > 4 )
1359 {
1360 // Pos phi diff
1361 phi1 = (*aRecMucTrackCol)[0]->getMucPos().phi(); phi2 = (*aRecMucTrackCol)[1]->getMucPos().phi();
1362 if( phi1 > 0 ) phiDiff = phi1 - phi2 - PI;
1363 else phiDiff = phi2 - phi1 - PI;
1364
1365 // Pos theta diff
1366 theta1 = (*aRecMucTrackCol)[0]->getMucPos().theta(); theta2 = (*aRecMucTrackCol)[1]->getMucPos().theta();
1367 thetaDiff = theta1 + theta2 - PI;
1368 m_hTrackPosPhiDiff->Fill( phiDiff );
1369 m_hTrackPosThetaDiff->Fill( thetaDiff );
1370 m_hDimuTracksPosDiff->Fill( thetaDiff, phiDiff );
1371 m_ntPosPhiDiff = phiDiff;
1372 m_ntPosThetaDiff = thetaDiff;
1373
1374 log << MSG::INFO << "PosPhiDiff:\t" << phiDiff << "\tPosThetaDiff:\t" << thetaDiff << endreq;
1375
1376 // Mom phi diff
1377 phi1 = (*aRecMucTrackCol)[0]->getMucMomentum().phi(); phi2 = (*aRecMucTrackCol)[1]->getMucMomentum().phi();
1378 if( phi1 > 0 ) phiDiff = phi1 - phi2 - PI;
1379 else phiDiff = phi2 - phi1 - PI;
1380
1381 // Mom theta diff
1382 theta1 = (*aRecMucTrackCol)[0]->getMucMomentum().theta(); theta2 = (*aRecMucTrackCol)[1]->getMucMomentum().theta();
1383 thetaDiff = theta1 + theta2 - PI;
1384
1385 m_hTrackMomPhiDiff->Fill( phiDiff );
1386 m_hTrackMomThetaDiff->Fill( thetaDiff );
1387 m_hDimuTracksMomDiff->Fill( thetaDiff, phiDiff );
1388 m_ntMomPhiDiff = phiDiff;
1389 m_ntMomThetaDiff = thetaDiff;
1390
1391 log << MSG::INFO << "MomPhiDiff:\t" << phiDiff << "\tMomThetaDiff:\t" << thetaDiff << endreq;
1392 m_ntDimuTag = m_eventTag;
1393 m_trackDiffTuple->write();
1394 }
1395
1396 // Retrieve tracks for calibration
1397 RecMucTrackCol::iterator trackIter = mucTrackCol->begin();
1398 int trackHitNum, rawHitNum, expectedHitNum, segNum, trkRecMode, lastLayerBR, lastLayerEC;
1399 int layerPassNum[3], passMax[TRACK_SEG_MAX][2];
1400 bool firedLay[TRACK_SEG_MAX][LAYER_MAX];
1401 bool seedList[PART_MAX][LAYER_MAX];
1402 trackHitNum = rawHitNum = expectedHitNum = segNum = trkRecMode = lastLayerBR = lastLayerEC = 0;
1403 layerPassNum[0] = layerPassNum[1] = layerPassNum[2] = 0;
1404 for( int segi=0; segi<TRACK_SEG_MAX; segi++ ) {
1405 passMax[segi][0] = passMax[segi][1] = 0;
1406 for( int layi=0; layi<LAYER_MAX; layi++ ) firedLay[segi][layi] = 0;
1407 }
1408
1409 mark_col trkSeg[TRACK_SEG_MAX];
1410
1411 // vector<MucRecHit*> mucRawHitCol;
1412 // vector<MucRecHit*> mucExpHitCol; //==========
1413 vector<int> mucRawHitCol;
1414 vector<int> mucExpHitCol;
1415
1416 for (int trackId = 0; trackIter != mucTrackCol->end(); trackIter++, trackId++)
1417 {
1418
1419//cout<<__LINE__<<"trackID "<<trackId<<endl;
1420 //trackHitNum = (*trackIter)->GetTotalHits();
1421 trackHitNum = (*trackIter)->numHits();
1422 log << MSG::DEBUG << "Track: " << trackId << " Hits: " << trackHitNum << endreq;
1423//cout<< "Track: " << trackId << " Hits: " << trackHitNum << endl;
1424 if( trackHitNum == 0 ) {
1425 log << MSG::INFO << "Track " << trackId << " no hits" << endreq;
1426 continue;
1427 }
1428
1429 m_ntTrackHits = trackHitNum;
1430
1431 m_trkRecMode = trkRecMode = (*trackIter)->GetRecMode();
1432 m_chi2 = (*trackIter)->chi2();
1433 m_px = (*trackIter)->getMucMomentum().x();
1434 m_py = (*trackIter)->getMucMomentum().y();
1435 m_pz = (*trackIter)->getMucMomentum().z();
1436 m_pt = sqrt(m_px*m_px + m_py*m_py);
1437 m_pp = sqrt(m_px*m_px + m_py*m_py + m_pz*m_pz);
1438
1439 // First fit position in MUC
1440 m_r = (*trackIter)->getMucPos().mag();
1441 m_cosTheta = (*trackIter)->getMucPos().cosTheta();
1442 m_theta = (*trackIter)->getMucPos().theta();
1443 m_phi = (*trackIter)->getMucPos().phi();
1444 m_depth = (*trackIter)->depth();
1445 m_brLastLayer = lastLayerBR = (*trackIter)->brLastLayer();
1446 m_ecLastLayer = lastLayerEC = (*trackIter)->ecLastLayer();
1447 m_totalHits = (*trackIter)->numHits();
1448 m_totalLayers = (*trackIter)->numLayers();
1449 m_maxHitsInLayer = (*trackIter)->maxHitsInLayer();
1450
1451
1452 m_hPhiCosTheta->Fill(m_cosTheta, m_phi);
1453 log << MSG::INFO << "Fill track info" << endreq;
1454
1455// MucRecHit* pMucRawHit;
1456// MucRecHit* pMucExpHit;
1457 if( m_calHitCol.size() != 0 ) m_calHitCol.clear(); // Fresh each track
1458
1459 // Digis belong to this rec track
1460 log << MSG::DEBUG << "Reconstruction hits(digis in a track): " << endreq;
1461 //mucRawHitCol = (*trackIter)->GetHits(); // Get hit collection of a track
1462 mucRawHitCol =(*trackIter)->getVecHits();
1463 rawHitNum += mucRawHitCol.size();
1464//cout<<__LINE__<<"rawHitNum "<<rawHitNum<<endl;
1465 segNum = 0;
1466/*
1467 if( trkRecMode == 3 ) { // By SlfTrk
1468 for(int iPart=0; iPart<PART_MAX; iPart++)
1469 for(int iLayer=0; iLayer<LAYER_MAX; iLayer++) seedList[iPart][iLayer] = false;
1470 }
1471*/
1472 for(unsigned int hitId = 0; hitId < mucRawHitCol.size(); hitId++)
1473 {
1474 /* pMucRawHit = mucRawHitCol[ hitId ];
1475 part = pMucRawHit->Part();
1476 segment = pMucRawHit->Seg();
1477 layer = pMucRawHit->Gap();
1478 strip = pMucRawHit->Strip();
1479 */
1480 mucId = mucRawHitCol[hitId];
1481 int part = MucID::barrel_ec(mucId);
1482 int segment = MucID::segment(mucId);
1483 int layer = MucID::layer(mucId);
1484 int strip = MucID::channel(mucId);
1485
1486//cout<<__LINE__<<" mu p s l s "<<mucId<<" "<<part<<" "<<segment<<" "<<layer<<" "<<strip<<endl;
1487 log << MSG::DEBUG << "[" << part << "\t" << segment << "\t" << layer << "\t" << strip << "]\t";
1488 //if( (hitId+1)%8 == 0 ) log << MSG::DEBUG << endreq;
1489
1490 // Add hit
1491 MucMark *aMark = new MucMark( part, segment, layer, strip );
1492 m_calHitCol.push_back( aMark );
1493
1494 // Set seed flag
1495// if( trkRecMode == 3 ) seedList[part][layer] = pMucRawHit->HitIsSeed(); // By SlfTrk
1496
1497 // Find track segment
1498 if(hitId == 0) { trkSeg[segNum].push_back( aMark ); segNum ++; }
1499 else
1500 {
1501 log << MSG::DEBUG << "segNum: " << segNum << endreq;
1502 bool notInSeg = true;
1503 for( int segi=0; segi<segNum; segi++ )
1504 {
1505 if( aMark->IsInSegWith( *(trkSeg[segi][0]) ) )
1506 {
1507 trkSeg[segi].push_back( aMark );
1508 notInSeg = false;
1509 break;
1510 }
1511 }
1512 // new track seg
1513 if( notInSeg == true )
1514 {
1515 trkSeg[segNum].push_back( aMark );
1516 segNum ++;
1517 if( segNum > TRACK_SEG_MAX ) {
1518 log << MSG::ERROR << "Track segment overflow: " << segNum << endreq;
1519 break;
1520 }
1521 }
1522 } // End else
1523 } // End raw hits
1524 log << MSG::DEBUG << endreq;
1525
1526 // Find maximal layers passed in track segments
1527 layerPassNum[0] = layerPassNum[1] = layerPassNum[2] = 0;
1528 for( int segi=0; segi<segNum; segi++ )
1529 {
1530 int tmpLayNum = 0;
1531 passMax[segi][0] = passMax[segi][1] = trkSeg[segi][0]->Layer();
1532 for( unsigned int hiti=1; hiti<trkSeg[segi].size(); hiti++ )
1533 {
1534 if( trkSeg[segi][hiti]->Layer() < passMax[segi][0] )
1535 passMax[segi][0] = trkSeg[segi][hiti]->Layer();
1536 if( trkSeg[segi][hiti]->Layer() > passMax[segi][1] )
1537 passMax[segi][1] = trkSeg[segi][hiti]->Layer();
1538 firedLay[segi][trkSeg[segi][hiti]->Layer()] = 1;
1539 }
1540
1541 for( int layi=0; layi<LAYER_MAX; layi++ ) {
1542 if( firedLay[segi][layi] ) tmpLayNum ++;
1543 }
1544
1545 if( segi == 0 ) layerPassNum[0] += passMax[segi][1] + 1;
1546 else layerPassNum[0] += (passMax[segi][1] - passMax[segi][0] + 1);
1547
1548 layerPassNum[1] += (passMax[segi][1] - passMax[segi][0] + 1);
1549 layerPassNum[2] += tmpLayNum;
1550
1551 trkSeg[segi].clear();
1552 }
1553 m_ntTrackEvent = m_currentEvent;
1554 m_ntTrackTag = m_eventTag;
1555 m_ntTrackSegFly = segNum;
1556 m_ntTrackLayFlyA = layerPassNum[0];
1557 m_ntTrackLayFlyB = layerPassNum[1];
1558 m_ntTrackLayFlyC = layerPassNum[2];
1559 m_trackInfoTuple->write();
1560 log << MSG::INFO << "Track\t" << trackId << "\tsegment(s):\t" << segNum
1561 << "\tlayer passed:\t" << layerPassNum[0] <<"\t" << layerPassNum[1] << "\t" << layerPassNum[2] << endreq;
1562 //if( layerPassNum[0]>B_LAY_NUM || layerPassNum[1]>B_LAY_NUM || layerPassNum[2]>B_LAY_NUM )
1563 // log << MSG::ERROR << "Over max layer:\t" << m_currentRun << "\t" << m_currentEvent << endreq;
1564
1565 // Expected hits in this rec track
1566 log << MSG::DEBUG << "Fitting hits(expected hits in a track): " << endreq;
1567 //mucExpHitCol = (*trackIter)->GetExpectedHits();
1568 mucExpHitCol = (*trackIter)->getExpHits();
1569
1570 expectedHitNum += mucExpHitCol.size();
1571//cout<<__LINE__<<"expectedHitNum "<<expectedHitNum<<endl;
1572 for(unsigned int hitId = 0; hitId < mucExpHitCol.size(); hitId++)
1573 {
1574// pMucRawHit = mucExpHitCol[ hitId ];
1575// part = pMucRawHit->Part(); segment = pMucRawHit->Seg();
1576// layer = pMucRawHit->Gap(); strip = pMucRawHit->Strip();
1577
1578 mucId =mucExpHitCol[hitId];
1579 int part = MucID::barrel_ec(mucId);
1580 int segment = MucID::segment(mucId);
1581 int layer = MucID::layer(mucId);
1582 int strip = MucID::channel(mucId);
1583//cout<<__LINE__<<part<<" "<<segment<<" "<<layer<<" "<<strip<<endl;
1584/*
1585 if( m_usePad != 0 )
1586 {
1587 pad = pMucRawHit->GetPadID(); padZ = pMucRawHit->GetIntersectZ();
1588 padX = pMucRawHit->GetIntersectY(); padY = pMucRawHit->GetIntersectX(); // Note: local coordinate
1589
1590 if( part != BRID )
1591 {
1592 if(segment == 1) { padX = -padX; }
1593 else if(segment == 2) { padX = -padX, padY = -padY; }
1594 else if(segment == 3) { padY = -padY; }
1595 }
1596 }
1597*/
1598 // Avoid bias in seed layers
1599 // if( seedList[part][layer] == true ) continue;
1600
1601 MucMark* currentMark = new MucMark( part, segment, layer, strip );
1602 m_expHitCol.push_back( currentMark );
1603 //log << MSG::DEBUG << "[" << part << "\t" << segment << "\t" << layer << "\t" << strip << "]\t";
1604 //if( (hitId+1)%8 == 0 ) log << MSG::DEBUG << endreq;
1605
1606 // Judge efficiency hit
1607 int isInPos = -1;
1608 bool isInEffWindow = false;
1609 isInPos = currentMark->IsInCol( m_segDigiCol[part][segment] );
1610
1611 // Avoid bias in outer layers caused by low momentum tracks
1612 if( part == BRID && (layer-lastLayerBR>1) ) continue;
1613 if( part != BRID && (layer-lastLayerEC>1) ) continue;
1614
1615 // Avoid bias in both sides of the innermost layer of Barrel
1616 if( part==BRID && layer==0 && (strip<2 || strip>45) )
1617 {
1618 if( isInPos != -1) // expHit is fired
1619 {
1620 m_record[part][segment][layer][strip][2] ++; // Efficiency hit number
1621 m_record[part][segment][layer][strip][1] ++; // Rec track number
1622 m_effHitCol.push_back( m_segDigiCol[part][segment][isInPos] );
1623
1624 if( m_usePad != 0 ) {
1625 m_h2DExpMap[part][segment][layer]->Fill(padX, padY);
1626 m_h2DHitMap[part][segment][layer]->Fill(padX, padY);
1627 }
1628
1629 }
1630 else {
1631 m_record[part][segment][layer][strip][1] ++;
1632 if( m_usePad != 0 ) m_h2DExpMap[part][segment][layer]->Fill(padX, padY);
1633 }
1634 continue; // Judge next hit
1635 }
1636
1637 // Eff calibration
1638 if( isInPos != -1 ) // expHit is fired
1639 {
1640 m_record[part][segment][layer][strip][2] ++; // Efficiency hit number
1641 m_record[part][segment][layer][strip][1] ++; // Rec track number
1642 m_effHitCol.push_back( m_segDigiCol[part][segment][isInPos] );
1643
1644 if( m_usePad != 0 ) {
1645 m_h2DExpMap[part][segment][layer]->Fill(padX, padY);
1646 m_h2DHitMap[part][segment][layer]->Fill(padX, padY);
1647 }
1648
1649 continue; // Judge next hit
1650 }
1651 else for(int tempStrip=0, hiti=-m_effWindow; hiti<=m_effWindow; hiti++ )
1652 {
1653 if( hiti == 0 ) continue;
1654 tempStrip = strip + hiti;
1655 if( tempStrip < 0 || tempStrip > m_ptrIdTr->GetStripMax(part,segment,layer) ) continue;
1656
1657 isInPos = m_ptrMucMark->IsInCol( part, segment, layer, tempStrip, m_segDigiCol[part][segment] );
1658 if( isInPos != -1 )
1659 {
1660 m_record[part][segment][layer][tempStrip][2] ++; // Efficiency hit number
1661 m_record[part][segment][layer][tempStrip][1] ++; // Rec track number
1662 m_effHitCol.push_back( m_segDigiCol[part][segment][isInPos] );
1663
1664 if( m_usePad != 0 ) {
1665 m_h2DExpMap[part][segment][layer]->Fill(padX, padY);
1666 m_h2DHitMap[part][segment][layer]->Fill(padX, padY);
1667 }
1668
1669 m_ntEffWindow = hiti;
1670 m_effWindowTuple->write();
1671 isInEffWindow = true;
1672 }
1673
1674 } // End else
1675
1676 if( isInEffWindow ) { continue; } // Judge next hit
1677 else { // A hit should be fired but not fired and not in the EffWindow
1678 m_record[part][segment][layer][strip][1] ++; // Rec track number
1679 if( m_usePad != 0 ) m_h2DExpMap[part][segment][layer]->Fill(padX, padY);
1680 }
1681
1682 } // End expected hits
1683
1684 // Fill residual, and for the other way of eff calculation
1685 log << MSG::INFO << "Fill residual" << endreq;
1686 vector<float> m_lineResCol = (*trackIter)->getDistHits();
1687 vector<float> m_quadResCol = (*trackIter)->getQuadDistHits();
1688 vector<float> m_extrResCol = (*trackIter)->getExtDistHits();
1689 int mode = (*trackIter)->GetRecMode();
1690
1691 for(unsigned int nres = 0; nres < m_lineResCol.size(); nres++ )
1692 if( fabs(m_lineResCol[nres])>resMax ) resMax = fabs(m_lineResCol[nres]);
1693
1694 log << MSG::INFO << "Good track for res" << endreq;
1695 if( trackHitNum > 4 && m_lineResCol[0] != -99) // track is good for res
1696 {
1697 // Fill res histograms
1698 bool firedFlag[PART_MAX][LAYER_MAX][2];
1699 for(int iprt=0; iprt<PART_MAX; iprt++)
1700 for(int jlay=0; jlay<LAYER_MAX; jlay++)
1701 firedFlag[iprt][jlay][0] = firedFlag[iprt][jlay][1] = false;
1702
1703 for(unsigned int hitId = 0; hitId < mucExpHitCol.size(); hitId++)
1704 {
1705// pMucExpHit = mucExpHitCol[ hitId ];
1706// part = pMucExpHit->Part(); segment = pMucExpHit->Seg(); layer = pMucExpHit->Gap();
1707
1708 mucId = mucExpHitCol[hitId];
1709
1710 int part = MucID::barrel_ec(mucId);
1711 int segment = MucID::segment(mucId);
1712 int layer = MucID::layer(mucId);
1713//cout<<__LINE__<<" m p s l"<<mucId<<" "<<part<<" "<<segment<<" "<<layer<<endl;
1714 firedFlag[part][layer][0] = true;
1715 }
1716
1717 log << MSG::INFO << "Fit res" << endreq;
1718 for(unsigned int hitId = 0; hitId < mucRawHitCol.size(); hitId++)
1719 {
1720 // pMucRawHit = mucRawHitCol[ hitId ];
1721 // part = pMucRawHit->Part(); segment = pMucRawHit->Seg(); layer = pMucRawHit->Gap();
1722
1723 mucId = mucRawHitCol[hitId];
1724 int part = MucID::barrel_ec(mucId);
1725 int segment = MucID::segment(mucId);
1726 int layer = MucID::layer(mucId);
1727//cout<<__LINE__<<" part "<< mucId<<" "<<part<<" "<<segment<<" "<<layer<<endl;
1728 if( part == BRID ) m_hBarrelResDist[layer]->Fill( m_lineResCol[hitId] );
1729 else m_hEndcapResDist[layer]->Fill( m_lineResCol[hitId] );
1730
1731 // if exp is true and fired is true, and not filled yet
1732 if( firedFlag[part][layer][0] == true && firedFlag[part][layer][1] == false )
1733 {
1734 m_resPart = part;
1735 m_resSegment = segment;
1736 m_resLayer = layer;
1737 m_lineRes = m_lineResCol[hitId];
1738 //m_quadRes = m_quadResCol[hitId];
1739 //m_extrRes = m_extrResCol[hitId];
1740 m_quadRes = 999;
1741 m_extrRes = 999;
1742 m_resFired = 1;
1743 m_resMode = mode;
1744 m_resInfoTuple->write();
1745//cout<<__LINE__<<" m_resInfoTuple->write "<<endl;
1746 }
1747
1748 firedFlag[part][layer][1] = true;
1749 }
1750
1751 log << MSG::INFO << "Exp res" << endreq;
1752 for(unsigned int hitId = 0; hitId < mucExpHitCol.size(); hitId++)
1753 {
1754// pMucExpHit = mucExpHitCol[ hitId ];
1755// part = pMucExpHit->Part(); segment = pMucExpHit->Seg(); layer = pMucExpHit->Gap();
1756
1757 mucId = mucExpHitCol[hitId];
1758 int part = MucID::barrel_ec(mucId);
1759 int segment = MucID::segment(mucId);
1760 int layer = MucID::layer(mucId);
1761//cout<<__LINE__<<mucId<<" "<<part<<" "<<segment<<" "<<layer<<endl;
1762 if(firedFlag[part][layer][0] == true && firedFlag[part][layer][1] == false)
1763 {
1764 m_resPart = part;
1765 m_resSegment = segment;
1766 m_resLayer = layer;
1767 m_lineRes = 1000;
1768 m_quadRes = 1000;
1769 m_extrRes = 1000;
1770 m_resFired = 0;
1771 m_resMode = mode;
1772 m_resInfoTuple->write();
1773 }
1774 }
1775
1776 } // End fill residual, if track is good for res
1777
1778 mucRawHitCol.clear();
1779 mucExpHitCol.clear();
1780
1781 } // End read all tracks
1782
1783 //if( resMax > 300 ) cout <<"Res too big!\t"<< m_fTotalEvent <<"\t"<< m_currentRun <<"\t"<< m_currentEvent <<"\t"<< resMax << endl;
1784
1785 m_ntTrackNum = mucTrackCol->size();
1786
1787 m_fTotalEffHit += rawHitNum;
1788 log << MSG::INFO << "Total hits in this event, raw: " << rawHitNum << "\texpected: " << expectedHitNum << endreq;
1789 //<--- End retrieve rec tracks
1790
1791 //---> Searching inc/noise hits
1792 log << MSG::INFO << "Searching inc/noise hits" << endreq;
1793 bool isNosHit;
1794 bool hasEffHit;
1795 for( unsigned int i=0; i < m_digiCol.size(); i++ )
1796 {
1797 isNosHit = true;
1798
1799 if( m_digiCol[i]->IsInCol( m_effHitCol ) !=-1) continue; // digi in effHitCol
1800 else
1801 {
1802 for( unsigned int j=0; j < m_clusterCol.size(); j++ )
1803 {
1804 hasEffHit = false;
1805 for( unsigned int k=0; k<m_clusterCol[j].size(); k++)
1806 {
1807 if( m_clusterCol[j][k]->IsInCol(m_effHitCol) != -1) // Clusters have efficiency hit
1808 {
1809 hasEffHit = true;
1810 break; // Out a cluster
1811 }
1812 }
1813
1814 if ( hasEffHit && (m_digiCol[i]->IsInCol( m_clusterCol[j] ) != -1) ) {
1815 isNosHit = false;
1816 break; // Out cluster collection
1817 }
1818 } // End cluster col
1819
1820 if( isNosHit ) {
1821 m_nosHitCol.push_back( m_digiCol[i] );
1822 m_fTotalNosHit ++;
1823 }
1824 }// End else
1825 } // End digi collection
1826
1827 return StatusCode::SUCCESS;
1828}
1829
1830// Check event
1831StatusCode MucCalibMgr::CheckEvent()
1832{
1833 MsgStream log(msgSvc, "MucCalibMgr");
1834 log << MSG::INFO << "Check event" << endreq;
1835
1836 // Find over digis in digis
1837 // Note that only one digi relates to one strip
1838 log << MSG::INFO << "Searching over digi" << endreq;
1839 int overDigi = 0;
1840 for( unsigned int i=0; i < m_digiCol.size(); i++ )
1841 for( unsigned int j=i+1; j < m_digiCol.size(); j++ ) {
1842 if( (*m_digiCol[i]) == (*m_digiCol[j]) ) overDigi ++;
1843 }
1844
1845 if( overDigi !=0 )
1846 log << MSG::ERROR << overDigi << " over digi found in digis!" << endreq;
1847
1848 // Find over hits in reconstruction hits
1849 // Note that two tracks may pass thought one strip
1850 log << MSG::INFO << "Searching over hit" << endreq;
1851 int overHit = 0;
1852 for( unsigned int i=0; i < m_expHitCol.size(); i++ )
1853 for( unsigned int j=i+1; j < m_expHitCol.size(); j++ ) {
1854 if( (*m_expHitCol[i]) == (*m_expHitCol[j]) ) overHit ++;
1855 }
1856
1857 if( overHit != 0 )
1858 log << MSG::WARNING << overHit << " over hits found in rec hits!" << endreq;
1859
1860 // Find hits of reconstruction not in MUCdigis
1861 log << MSG::INFO << "Searching new hit" << endreq;
1862 int newHit = 0;
1863 int num = 0;
1864 for( unsigned int i=0; i < m_expHitCol.size(); i++ )
1865 {
1866 num = m_expHitCol[i]->NumInCol( m_digiCol );
1867 if( num == 0 )
1868 {
1869 log << MSG::ERROR << "Exp hit not in digis!"
1870 << "prt: " << (*m_expHitCol[i]).Part()
1871 << "\tseg: " << (*m_expHitCol[i]).Segment()
1872 << "\tlay: " << (*m_expHitCol[i]).Layer()
1873 << "\tstr: " << (*m_expHitCol[i]).Strip() << endreq;
1874
1875 newHit ++;
1876 }
1877 }
1878
1879 if( newHit != 0 )
1880 log << MSG::WARNING << newHit << " new hit(s) found in rec hits!" << endreq;
1881
1882 return StatusCode::SUCCESS;
1883}
1884
1885// Fill event
1886StatusCode MucCalibMgr::FillEvent()
1887{
1888 MsgStream log(msgSvc, "MucCalibMgr");
1889 log << MSG::INFO << "Fill event" << endreq;
1890
1891 int part, segment, layer, strip, size;
1892 part = segment = layer = strip = size = 0;
1893
1894 // Fill digis
1895 log << MSG::INFO << "Fill digis" << endreq;
1896 for( unsigned int i=0; i<m_digiCol.size(); i++ )
1897 {
1898 part = (*m_digiCol[i]).Part();
1899 segment = (*m_digiCol[i]).Segment();
1900 layer = (*m_digiCol[i]).Layer();
1901 strip = (*m_digiCol[i]).Strip();
1902
1903 FillDigi( part, segment, layer, strip );
1904 m_record[part][segment][layer][strip][0] ++;
1905 m_fTotalDigi ++;
1906 }
1907
1908 // Fill rec hits
1909 log << MSG::INFO << "Fill rec hits" << endreq;
1910 for( unsigned int i=0; i<m_expHitCol.size(); i++ )
1911 {
1912 part = (*m_expHitCol[i]).Part();
1913 segment = (*m_expHitCol[i]).Segment();
1914 layer = (*m_expHitCol[i]).Layer();
1915 strip = (*m_expHitCol[i]).Strip();
1916
1917 FillExpHit( part, segment, layer, strip );
1918 m_record[part][segment][layer][strip][4] ++;
1919 m_fTotalExpHit ++;
1920 }
1921
1922 // Fill eff hits
1923 log << MSG::INFO << "Fill eff hits" << endreq;
1924 for( unsigned int i=0; i<m_effHitCol.size(); i++ )
1925 {
1926 part = (*m_effHitCol[i]).Part();
1927 segment = (*m_effHitCol[i]).Segment();
1928 layer = (*m_effHitCol[i]).Layer();
1929 strip = (*m_effHitCol[i]).Strip();
1930
1931 FillEffHit( part, segment, layer, strip );
1932 m_fTotalEffHit ++;
1933 }
1934
1935 // Fill clusters
1936 log << MSG::INFO << "Fill clusters" << endreq;
1937 for( unsigned int i=0; i<m_clusterCol.size(); i++ )
1938 {
1939 size = m_clusterCol[i].size();
1940 // may include the special cluster, size = 1
1941 if( size > CLUSTER_CUT )
1942 {
1943 part = (*m_clusterCol[i][0]).Part();
1944 segment = (*m_clusterCol[i][0]).Segment();
1945 layer = (*m_clusterCol[i][0]).Layer();
1946
1947 FillCluster( part, segment, layer, size );
1948 m_ntClusterSize = size;
1949 m_clusterSizeTuple->write();
1950 }
1951 }
1952
1953 // Fill inc/noise hits
1954 log << MSG::INFO << "Fill inc/noise hits" << endreq;
1955 for( unsigned int i=0; i<m_nosHitCol.size(); i++ )
1956 {
1957 part = (*m_nosHitCol[i]).Part();
1958 segment = (*m_nosHitCol[i]).Segment();
1959 layer = (*m_nosHitCol[i]).Layer();
1960 strip = (*m_nosHitCol[i]).Strip();
1961
1962 FillNosHit( part, segment, layer, strip );
1963 m_record[part][segment][layer][strip][3] ++;
1964 }
1965
1966 // Event log
1967 m_ntEventId = m_currentEvent;
1968 m_ntEventTag = m_eventTag;
1969 m_ntDigiNum = m_digiCol.size();
1970 m_ntExpHitNum = m_expHitCol.size();
1971 m_ntEffHitNum = m_effHitCol.size();
1972 m_ntNosHitNum = m_nosHitCol.size();
1973 m_ntClusterNum = m_clusterCol.size();
1974
1975 // Reset mark collection
1976 for( unsigned int i=0; i<m_digiCol.size(); i++ ) {
1977 if( m_digiCol[i] != NULL ) delete m_digiCol[i];
1978 }
1979
1980 for( unsigned int i=0; i<m_expHitCol.size(); i++ ) {
1981 if( m_expHitCol[i] != NULL ) delete m_expHitCol[i];
1982 }
1983/*
1984 for( unsigned int i=0; i<m_effHitCol.size(); i++ ) {
1985 if( m_effHitCol[i] != NULL ) delete m_effHitCol[i];
1986 }
1987 log << MSG::INFO << m_effHitCol.size() << endreq;
1988*/
1989 for( unsigned int i=0; i<m_calHitCol.size(); i++ ) {
1990 if( m_calHitCol[i] != NULL ) delete m_calHitCol[i];
1991 }
1992
1993 if( m_digiCol.size() != 0 ) m_digiCol.clear();
1994 if( m_expHitCol.size() != 0 ) m_expHitCol.clear();
1995 if( m_calHitCol.size() != 0 ) m_calHitCol.clear();
1996 if( m_effHitCol.size() != 0 ) m_effHitCol.clear();
1997 if( m_nosHitCol.size() != 0 ) m_nosHitCol.clear();
1998 if( m_clusterCol.size() != 0 ) m_clusterCol.clear();
1999
2000
2001 for( int i=0; i<PART_MAX; i++ ) {
2002 for( int j=0; j<SEGMENT_MAX; j++ ) {
2003 if( m_segDigiCol[i][j].size() != 0 ) m_segDigiCol[i][j].clear();
2004 }
2005 }
2006
2007 m_evtEnd = clock();
2008
2009 m_ntEventTime = (double(m_evtEnd - m_evtBegin))/CLOCKS_PER_SEC;
2010 log << MSG::INFO << "Event time:\t" << m_ntEventTime << "s" << endreq;
2011 m_eventLogTuple->write();
2012
2013 return StatusCode::SUCCESS;
2014}
2015
2016StatusCode MucCalibMgr::FillDigi( int part, int segment, int layer, int strip )
2017{
2018 // Hit map for online check
2019 int tmpId;
2020
2021 if( (int)m_tag || m_dimuOnly==0 || (m_dimuOnly==1&&m_eventTag==1) )
2022 {
2023 if( part == BRID )
2024 {
2025 // According to layer
2026 if( layer%2 == 0 ) {
2027 if( segment > 4 ) tmpId = segment*48 + (47 - strip);
2028 else tmpId = segment*48 + strip;
2029 }
2030 else if( segment < 3 ) tmpId = segment*96 + strip;
2031 else tmpId = (segment-1)*96 + 112 + strip;
2032
2033 m_hHitMapBarrel_Lay[layer]->Fill(tmpId);
2034
2035 // According to segment
2036 if( segment != B_TOP ) {
2037 if( segment > 4 )
2038 tmpId = 48*((layer+1)/2) + 96*(layer/2) + ( ((layer+1)%2)*48 + (layer%2)*96 - strip - 1 );
2039 else tmpId = 48*((layer+1)/2) + 96*(layer/2) + strip;
2040 }
2041 else tmpId = 48*((layer+1)/2) + 112*(layer/2) + strip;
2042
2043 m_hHitMapBarrel_Seg[segment]->Fill(tmpId);
2044 }
2045 else if( part == EEID )
2046 {
2047 // According to layer
2048 tmpId = segment*64 + strip;
2049 m_hHitMapEndcap_Lay[0][layer]->Fill(tmpId);
2050 // According to segment
2051 tmpId = layer*64 + strip;
2052 m_hHitMapEndcap_Seg[0][segment]->Fill(tmpId);
2053 }
2054 else if( part == EWID )
2055 {
2056 // According to layer
2057 tmpId = segment*64 + strip;
2058 m_hHitMapEndcap_Lay[1][layer]->Fill(tmpId);
2059 // According to segment
2060 tmpId = layer*64 + strip;
2061 m_hHitMapEndcap_Seg[1][segment]->Fill(tmpId);
2062 }
2063 }
2064
2065 // Total units
2066 int boxId = m_ptrIdTr->GetBoxId( part, segment, layer );
2067 int strId = m_ptrIdTr->GetStripId( part, segment, layer, strip );
2068
2069 m_hStripFireMap[boxId]->Fill( strip );
2070 m_hStripFire->Fill( strId );
2071 m_hBoxFire->Fill( boxId );
2072 m_hLayerFire->Fill( layer );
2073 if( part==BRID ) m_hBrLayerFire->Fill( layer );
2074 else if( part== EEID ) m_hEcLayerFire->Fill( layer+1 );
2075 else m_hEcLayerFire->Fill( -(layer+1) );
2076
2077 return StatusCode::SUCCESS;
2078}
2079
2080StatusCode MucCalibMgr::FillExpHit( int part, int segment, int layer, int strip )
2081{
2082 int boxId = m_ptrIdTr->GetBoxId( part, segment, layer );
2083 int strId = m_ptrIdTr->GetStripId( part, segment, layer, strip );
2084
2085 m_hStripExpHitMap[boxId]->Fill( strip );
2086 m_hStripExpHit->Fill( strId );
2087 m_hBoxExpHit->Fill( boxId );
2088 m_hLayerExpHit->Fill( layer );
2089
2090 return StatusCode::SUCCESS;
2091}
2092
2093StatusCode MucCalibMgr::FillEffHit( int part, int segment, int layer, int strip )
2094{
2095 int boxId = m_ptrIdTr->GetBoxId( part, segment, layer );
2096 int strId = m_ptrIdTr->GetStripId( part, segment, layer, strip );
2097
2098 m_hStripEffHitMap[boxId]->Fill( strip );
2099 m_hStripEffHit->Fill( strId );
2100 m_hBoxEffHit->Fill( boxId );
2101 m_hLayerEffHit->Fill( layer );
2102
2103 return StatusCode::SUCCESS;
2104}
2105
2106StatusCode MucCalibMgr::FillNosHit( int part, int segment, int layer, int strip )
2107{
2108 int boxId = m_ptrIdTr->GetBoxId( part, segment, layer );
2109 int strId = m_ptrIdTr->GetStripId( part, segment, layer, strip );
2110
2111 m_hStripNosHitMap[boxId]->Fill( strip );
2112 m_hStripNosHit->Fill( strId );
2113 m_hBoxNosHit->Fill( boxId );
2114 m_hLayerNosHit->Fill( layer );
2115
2116 return StatusCode::SUCCESS;
2117}
2118
2119StatusCode MucCalibMgr::FillCluster( int part, int segment, int layer, int size )
2120{
2121 int boxId = m_ptrIdTr->GetBoxId( part, segment, layer );
2122
2123 m_hLayerCluster[layer]->Fill( size );
2124 m_hBoxCluster[boxId]->Fill( size );
2125
2126 return StatusCode::SUCCESS;
2127}
2128
2129
2130//-------------------------------------------------------------------------------------------------
2131//-------------------------------- Member functions for finalizing --------------------------------
2132//-------------------------------------------------------------------------------------------------
2133
2134StatusCode MucCalibMgr::AnalyseEffAndNoise()
2135{
2136 MsgStream log(msgSvc, "MucCalibMgr");
2137 log<< MSG::INFO << "Start efficiency and noise calibration" << endreq;
2138
2139 // Set time
2140 m_fTotalDAQTime = m_fTotalEvent/m_er;
2141 log<< MSG::INFO << "Total run time:\t" << m_fTotalDAQTime << endreq;
2142
2143 EffAndNoiseLV2();
2144 EffAndNoiseLV1();
2145 EffAndNoiseLV0();
2146
2147 if( m_usePad != 0 ) PadEff();
2148
2149 return StatusCode::SUCCESS;
2150}
2151
2152StatusCode MucCalibMgr::EffAndNoiseLV0()
2153{
2154 MsgStream log(msgSvc, "MucCalibMgr");
2155 log<< MSG::INFO << "Analyse layer efficiency and noise" << endreq;
2156
2157 int part, segment, layer, stripMax;
2158 part = segment = layer = stripMax = 0;
2159
2160 double digi, effHit, trackNum, nosHit, recHit;
2161 double eff, effErr, noise, nosRatio, nosRatioErr, cnt, cluster;
2162 eff = effErr = noise = nosRatio = nosRatioErr = cnt = cluster = 0.;
2163
2164 for( int i=0; i<LAYER_MAX; i++ )
2165 {
2166 digi = effHit = trackNum = nosHit = recHit = 0;
2167
2168 for( int j=0; j<PART_MAX; j++ )
2169 {
2170 for( int k=0; k<SEGMENT_MAX; k++)
2171 {
2172 stripMax = m_ptrIdTr->GetStripMax( j, k, i );
2173 for( int n=0; n<stripMax; n++ )
2174 {
2175 digi += m_record[j][k][i][n][0];
2176 trackNum += m_record[j][k][i][n][1];
2177 effHit += m_record[j][k][i][n][2];
2178 nosHit += m_record[j][k][i][n][3];
2179 recHit += m_record[j][k][i][n][4];
2180 }
2181 }
2182 }
2183
2184 // Efficiency
2185 if( trackNum == 0 ) {
2186 eff = effErr = 0.0;
2187 //eff = DEFAULT_EFF_VALUE;
2188 //effErr = DEFAULT_EFF_ERR;
2189 }
2190 else
2191 {
2192 eff = ( (double)effHit )/trackNum;
2193 effErr = sqrt( eff*(1-eff)/trackNum );
2194 m_fCalibLayerNum ++;
2195 }
2196
2197 // Noise
2198 if( m_layerResults[3][i] < LIMIT_CUT || m_fTotalDAQTime < LIMIT_CUT )
2199 noise = DEFAULT_NOS_VALUE;
2200 else {
2201 if( m_recMode == 2 ) noise = (double)nosHit/(m_fTotalEvent*TRIGGER_WINDOW*m_layerResults[3][i]);
2202 else noise = (double)nosHit/(m_fTotalDAQTime * m_layerResults[3][i]);
2203 }
2204
2205 if( digi != 0 ) {
2206 nosRatio = ( (double)nosHit )/digi;
2207 nosRatioErr = sqrt( nosRatio*(1-nosRatio)/digi );
2208 }
2209 else {
2210 nosRatio = DEFAULT_INC_VALUE;
2211 nosRatioErr = 0;
2212 }
2213
2214 // Counting rate
2215 if( m_recMode == 2 )
2216 cnt = (double)digi/(m_fTotalEvent*TRIGGER_WINDOW * m_layerResults[3][i]);
2217 else
2218 cnt = (double)digi/(m_fTotalDAQTime * m_layerResults[3][i]);
2219
2220 // Cluster
2221 cluster = m_hLayerCluster[i]->GetMean();
2222
2223 m_layerResults[0][ i ] = eff;
2224 m_layerResults[1][ i ] = effErr;
2225 m_layerResults[2][ i ] = noise;
2226 m_layerResults[4][ i ] = cluster;
2227 m_layerResults[5][ i ] = trackNum;
2228
2229 // Fill histogram
2230 m_hLayerEff->Fill( i, eff );
2231 m_hLayerEff->SetBinError( i+1, effErr );
2232 m_hLayerNosRatio->Fill( i, nosRatio );
2233 m_hLayerNosRatio->SetBinError( i+1, nosRatioErr );
2234 m_hLayerNos->Fill( i, noise );
2235 m_hLayerCnt->Fill( i, cnt );
2236
2237 // Add cluster histogram
2238 m_hLayerClusterCmp->Fill( i, cluster );
2239
2240 // Fill tree
2241 m_fLayerId = i;
2242 m_fLayerEff = eff;
2243 m_fLayerEffErr = effErr;
2244 m_fLayerTrkNum = trackNum;
2245 m_fLayerExpHit = recHit;
2246 m_fLayerEffHit = effHit;
2247 m_fLayerNosRatio = nosRatio;
2248 m_fLayerDigi = digi;
2249 m_fLayerNosHit = nosHit;
2250 m_fLayerNos = noise;
2251 m_fLayerCnt = cnt;
2252 m_tLayConst->Fill();
2253
2254 // Add cluster ntuple
2255 m_fLayerCluster = cluster;
2256
2257 } // End layer
2258
2259 return StatusCode::SUCCESS;
2260}
2261
2262StatusCode MucCalibMgr::EffAndNoiseLV1()
2263{
2264 MsgStream log(msgSvc, "MucCalibMgr");
2265 log<< MSG::INFO << "Analyse box efficiency and noise" << endreq;
2266
2267 int part, segment, layer, stripMax;
2268 part = segment = layer = stripMax = 0;
2269
2270 double digi, effHit, trackNum, nosHit, recHit;
2271 double eff, effErr, noise, nosRatio, nosRatioErr, cnt, cluster;
2272 eff = effErr = noise = nosRatio = nosRatioErr = cnt = cluster = 0.;
2273
2274 for( int i=0; i<BOX_MAX; i++ )
2275 {
2276 m_ptrIdTr->SetBoxPos( i, &part, &segment, &layer );
2277 stripMax = m_ptrIdTr->GetStripMax( part, segment, layer );
2278
2279 digi = effHit = trackNum = nosHit = recHit = 0;
2280 for( int j=0; j<stripMax; j++ )
2281 {
2282 digi += m_record[part][segment][layer][j][0];
2283 trackNum += m_record[part][segment][layer][j][1];
2284 effHit += m_record[part][segment][layer][j][2];
2285 nosHit += m_record[part][segment][layer][j][3];
2286 recHit += m_record[part][segment][layer][j][4];
2287 }
2288
2289 // Efficiency
2290 if( trackNum == 0 ) {
2291 eff = effErr = 0.0;
2292 //eff = DEFAULT_EFF_VALUE;
2293 //effErr = DEFAULT_EFF_ERR;
2294 }
2295 else
2296 {
2297 eff = ( (double)effHit )/trackNum;
2298 effErr = sqrt( eff*(1-eff)/trackNum );
2299 m_fCalibBoxNum ++;
2300 }
2301
2302 // Noise
2303 if( m_boxResults[3][i] < LIMIT_CUT || m_fTotalDAQTime < LIMIT_CUT )
2304 noise = DEFAULT_NOS_VALUE;
2305 else {
2306 if( m_recMode == 2 )
2307 noise = (double)nosHit/(m_fTotalEvent*TRIGGER_WINDOW * m_boxResults[3][i]);
2308 else
2309 noise = (double)nosHit/(m_fTotalDAQTime * m_boxResults[3][i]);
2310 }
2311
2312 if( digi != 0 ) {
2313 nosRatio = ( (double)nosHit )/digi;
2314 nosRatioErr = sqrt( nosRatio*(1-nosRatio)/digi );
2315 }
2316 else {
2317 nosRatio = DEFAULT_INC_VALUE;
2318 nosRatioErr = 0;
2319 }
2320
2321 // Counting rate
2322 if( m_recMode == 2 )
2323 cnt = (double)digi/(m_fTotalEvent*TRIGGER_WINDOW * m_boxResults[3][i]);
2324 else
2325 cnt = (double)digi/(m_fTotalDAQTime * m_boxResults[3][i]);
2326
2327 // Cluster
2328 cluster = m_hBoxCluster[i]->GetMean();
2329
2330 m_boxResults[0][ i ] = eff;
2331 m_boxResults[1][ i ] = effErr;
2332 m_boxResults[2][ i ] = noise;
2333 m_boxResults[4][ i ] = cluster;
2334 m_boxResults[5][ i ] = trackNum;
2335
2336 // Fill histograms
2337 m_hBoxEff->Fill( i, eff );
2338 m_hBoxEff->SetBinError( i+1, effErr );
2339 m_hBoxNosRatio->Fill( i, nosRatio );
2340 m_hBoxNosRatio->SetBinError( i+1, nosRatioErr );
2341 m_hBoxNos->Fill( i, noise );
2342 m_hBoxCnt->Fill( i, cnt );
2343
2344 // add cluster histogram
2345 m_hBoxClusterCmp->Fill( i, cluster );
2346
2347 // Fill tree
2348 m_fBoxId = i;
2349 m_fBoxPart = part;
2350 m_fBoxSegment = segment;
2351 m_fBoxLayer = layer;
2352 m_fBoxEff = eff;
2353 m_fBoxEffErr = effErr;
2354 m_fBoxTrkNum = trackNum;
2355 m_fBoxExpHit = recHit;
2356 m_fBoxEffHit = effHit;
2357 m_fBoxNosRatio = nosRatio;
2358 m_fBoxDigi = digi;
2359 m_fBoxNosHit = nosHit;
2360 m_fBoxNos = noise;
2361 m_fBoxCnt = cnt;
2362 m_tBoxConst->Fill();
2363
2364 // add box cluster ntuple
2365 m_fBoxCluster = cluster;
2366
2367 }// End BOX_MAX
2368
2369 return StatusCode::SUCCESS;
2370}
2371
2372StatusCode MucCalibMgr::EffAndNoiseLV2()
2373{
2374 MsgStream log(msgSvc, "MucCalibMgr");
2375 log<< MSG::INFO << "Analyse strip efficiency and noise" << endreq;
2376
2377 int part, segment, layer, stripMax;
2378 part = segment = layer = stripMax = 0;
2379
2380 double digi, effHit, trackNum, nosHit, recHit;
2381 double eff, effErr, noise, nosRatio, nosRatioErr, cnt, cluster;
2382 eff = effErr = noise = nosRatio = nosRatioErr = cnt = cluster = 0.;
2383
2384 for( int i=0; i<BOX_MAX; i++ )
2385 {
2386 m_ptrIdTr->SetBoxPos( i, &part, &segment, &layer );
2387 stripMax = m_ptrIdTr->GetStripMax( part, segment, layer );
2388
2389 for( int j=0; j<stripMax; j++ )
2390 {
2391 digi = m_record[part][segment][layer][j][0];
2392 trackNum = m_record[part][segment][layer][j][1];
2393 effHit = m_record[part][segment][layer][j][2];
2394 nosHit = m_record[part][segment][layer][j][3];
2395 recHit = m_record[part][segment][layer][j][4];
2396
2397 int stripId = m_ptrIdTr->GetStripId( part, segment, layer, j );
2398
2399 // Efficiency
2400 if( trackNum == 0 ) {
2401 eff = effErr = 0.0;
2402 //eff = DEFAULT_EFF_VALUE;
2403 //effErr = DEFAULT_EFF_ERR;
2404 }
2405 else
2406 {
2407 eff = ( (double)effHit )/trackNum;
2408 effErr = sqrt( eff*(1-eff)/trackNum );
2409 m_fCalibStripNum ++;
2410 }
2411
2412 // Noise
2413 if( m_stripResults[3][stripId] < LIMIT_CUT || m_fTotalDAQTime < LIMIT_CUT )
2414 noise = DEFAULT_NOS_VALUE;
2415 else {
2416 if( m_recMode == 2 )
2417 noise = (double)nosHit/(m_fTotalEvent*TRIGGER_WINDOW * m_stripResults[3][stripId]);
2418 else
2419 noise = (double)nosHit/(m_fTotalDAQTime * m_stripResults[3][stripId]);
2420 }
2421
2422 if( digi != 0 ) {
2423 nosRatio = ( (double)nosHit )/digi;
2424 nosRatioErr = sqrt( nosRatio*(1-nosRatio)/digi );
2425 }
2426 else {
2427 nosRatio = DEFAULT_INC_VALUE;
2428 nosRatioErr = 0.;
2429 }
2430
2431 // Counting rate
2432 if( m_recMode == 2 )
2433 cnt = (double)digi/(m_fTotalEvent*TRIGGER_WINDOW * m_stripResults[3][stripId]);
2434 else
2435 cnt = (double)digi/(m_fTotalDAQTime * m_stripResults[3][stripId]);
2436
2437 // Strip itself no clusters
2438 m_stripResults[0][ stripId ] = eff;
2439 m_stripResults[1][ stripId ] = effErr;
2440 m_stripResults[2][ stripId ] = noise;
2441 m_stripResults[5][ stripId ] = trackNum;
2442
2443 // Fill histotram
2444 m_hStripEffMap[i]->Fill( j, eff );
2445 m_hStripEffMap[i]->SetBinError( j+1, effErr );
2446 m_hStripEff->Fill( stripId, eff );
2447 m_hStripEff->SetBinError( stripId+1, effErr );
2448 m_hStripNosRatioMap[i]->Fill( j, nosRatio );
2449 m_hStripNosRatioMap[i]->SetBinError( j+1, nosRatioErr );
2450 m_hStripNosRatio->Fill( stripId, nosRatio );
2451 m_hStripNosRatio->SetBinError( stripId+1, nosRatioErr );
2452 m_hStripNos->Fill( stripId, noise );
2453 m_hStripCnt->Fill( stripId, cnt );
2454
2455 // Fill Tree
2456 m_fStripId = stripId;
2457 m_fStripPart = part;
2458 m_fStripSegment = segment;
2459 m_fStripLayer = layer;
2460 m_fStripEff = eff;
2461 m_fStripEffErr = effErr;
2462 m_fStripNosRatio = nosRatio;
2463 m_fStripDigi = digi;
2464 m_fStripNos = noise;
2465 m_fStripCnt = cnt;
2466 m_fStripEffHit = effHit;
2467 m_fStripExpHit = recHit;
2468 m_fStripNosHit = nosHit;
2469 m_fStripTrkNum = trackNum;
2470 m_tStrConst->Fill();
2471
2472 } // End stripMax
2473 } // End BOX_MAX
2474
2475 return StatusCode::SUCCESS;
2476}
2477
2478StatusCode MucCalibMgr::PadEff()
2479{
2480 MsgStream log(msgSvc, "MucCalibMgr");
2481
2482 int xBinStart, xBinEnd, yBinStart, yBinEnd;
2483 double expHit, firedHit, eff;
2484 eff = expHit = firedHit = 0.;
2485
2486 for( int i=0; i<PART_MAX; i++ )
2487 for( int j=0; j<((i==BRID)?B_SEG_NUM:E_SEG_NUM); j++ )
2488 for( int k=0; k<((i==BRID)?B_LAY_NUM:E_LAY_NUM); k++ )
2489 {
2490 xBinStart = m_h2DExpMap[i][j][k]->GetXaxis()->GetFirst();
2491 xBinEnd = m_h2DExpMap[i][j][k]->GetXaxis()->GetLast() + 1;
2492 yBinStart = m_h2DExpMap[i][j][k]->GetYaxis()->GetFirst();
2493 yBinEnd = m_h2DExpMap[i][j][k]->GetYaxis()->GetLast() + 1;
2494
2495 for( int xi = xBinStart; xi<xBinEnd; xi++ )
2496 for( int yi = yBinStart; yi<yBinEnd; yi++ )
2497 {
2498 expHit = m_h2DExpMap[i][j][k]->GetBinContent(xi, yi);
2499 firedHit = m_h2DHitMap[i][j][k]->GetBinContent(xi, yi);
2500
2501 if( expHit !=0 ) eff = firedHit / expHit;
2502 else eff = 0;
2503
2504 if( eff>1.0 )
2505 cout<<"Eff error:\t["<<i<<"\t"<<j<<"\t"<<k<<",\t"<<xi<<"\t"<<yi<<"]:\t"
2506 <<eff<<"\t,"<<firedHit<<" "<<expHit<<endl;
2507
2508 m_h2DEffMap[i][j][k]->SetBinContent(xi, yi, eff);
2509 }
2510 }
2511 return StatusCode::SUCCESS;
2512}
2513
2514StatusCode MucCalibMgr::AnalyseCluster()
2515{
2516 MsgStream log(msgSvc, "MucCalibMgr");
2517
2518 return StatusCode::SUCCESS;
2519}
2520
2521StatusCode MucCalibMgr::AnalyseRes()
2522{
2523 MsgStream log(msgSvc, "MucCalibMgr");
2524 log << MSG::INFO << "Analyse spacial resolution" << endreq;
2525
2526 double resMean, resMeanErr, resSigma, resSigmaErr;
2527 resMean = resMeanErr = resSigma = resSigmaErr = 0.;
2528
2529 for( int i=0; i<B_LAY_NUM; i++ )
2530 {
2531 m_hBarrelResDist[i]->Fit("gaus");
2532 resMean = m_hBarrelResDist[i]->GetFunction("gaus")->GetParameter("Mean");
2533 resSigma = m_hBarrelResDist[i]->GetFunction("gaus")->GetParameter("Sigma");
2534 resMeanErr = m_hBarrelResDist[i]->GetFunction("gaus")->GetParError(1);
2535 resSigmaErr = m_hBarrelResDist[i]->GetFunction("gaus")->GetParError(2);
2536
2537 m_hBarrelResComp[0]->Fill( i, resMean );
2538 m_hBarrelResComp[0]->SetBinError( i+1, resMeanErr );
2539 m_hBarrelResComp[1]->Fill( i, resSigma );
2540 m_hBarrelResComp[1]->SetBinError( i+1, resSigmaErr );
2541 }
2542
2543 for( int i=0; i<E_LAY_NUM; i++ )
2544 {
2545 m_hEndcapResDist[i]->Fit("gaus");
2546 resMean = m_hEndcapResDist[i]->GetFunction("gaus")->GetParameter("Mean");
2547 resSigma = m_hEndcapResDist[i]->GetFunction("gaus")->GetParameter("Sigma");
2548 resMeanErr = m_hEndcapResDist[i]->GetFunction("gaus")->GetParError(1);
2549 resSigmaErr = m_hEndcapResDist[i]->GetFunction("gaus")->GetParError(2);
2550
2551 m_hEndcapResComp[0]->Fill( i, resMean );
2552 m_hEndcapResComp[0]->SetBinError( i+1, resMeanErr );
2553 m_hEndcapResComp[1]->Fill( i, resSigma );
2554 m_hEndcapResComp[1]->SetBinError( i+1, resSigmaErr );
2555 }
2556
2557 return StatusCode::SUCCESS;
2558}
2559
2560StatusCode MucCalibMgr::SaveConst()
2561{
2562 MsgStream log(msgSvc, "MucCalibMgr");
2563 log << MSG::INFO << "Save calibration constants" << endreq;
2564
2565 // Save calibrated eff in graphes
2566 // LV0
2567 double layerXY[2][LAYER_MAX];
2568 double layerEXY[2][LAYER_MAX];
2569 for( int i=0; i<LAYER_MAX; i++ )
2570 {
2571 layerXY[0][i] = i;
2572 layerEXY[0][i] = 0;
2573 if( m_layerResults[5][i] >= 100*TRACK_THRESHOLD ) {
2574 layerXY[1][i] = m_layerResults[0][i];
2575 layerEXY[1][i] = m_layerResults[1][i];
2576 }
2577 else {
2578 layerXY[1][i] = 0;
2579 layerEXY[1][i] = 0;
2580 }
2581 }
2582 m_geLayerEff = new TGraphErrors(LAYER_MAX, layerXY[0], layerXY[1], layerEXY[0], layerEXY[1]);
2583 m_geLayerEff->SetMarkerStyle(25);
2584 m_geLayerEff->SetMarkerSize(0.5);
2585 m_cv[0] = new TCanvas("GoodLayerEff","Layer efficiency", 50, 50, 800, 600);
2586 m_cv[0]->SetFillColor(0);
2587 m_cv[0]->SetBorderMode(0);
2588 m_geLayerEff->Draw("AP");
2589 m_cv[0]->Write();
2590
2591 // LV1
2592 double boxXY[2][BOX_MAX];
2593 double boxEXY[2][BOX_MAX];
2594 for( int i=0; i<BOX_MAX; i++ )
2595 {
2596 boxXY[0][i] = i;
2597 boxEXY[0][i] = 0;
2598 if( m_boxResults[5][i] >= 10*TRACK_THRESHOLD ) {
2599 boxXY[1][i] = m_boxResults[0][i];
2600 boxEXY[1][i] = m_boxResults[1][i];
2601 }
2602 else {
2603 boxXY[1][i] = 0;
2604 boxEXY[1][i] = 0;
2605 }
2606 }
2607 m_geBoxEff = new TGraphErrors(BOX_MAX, boxXY[0], boxXY[1], boxEXY[0], boxEXY[1]);
2608 m_geBoxEff->SetMarkerStyle(25);
2609 m_geBoxEff->SetMarkerSize(0.5);
2610 m_cv[1] = new TCanvas("GoodBoxEff","Box efficiency", 75, 75, 800, 600);
2611 m_cv[1]->SetFillColor(0);
2612 m_cv[1]->SetBorderMode(0);
2613 m_geBoxEff->Draw("AP");
2614 m_cv[1]->Write();
2615
2616 // LV2
2617 double stripXY[2][STRIP_MAX];
2618 double stripEXY[2][STRIP_MAX];
2619 for( int i=0; i<STRIP_MAX; i++ )
2620 {
2621 stripXY[0][i] = i;
2622 stripEXY[0][i] = 0;
2623 if( m_stripResults[5][i] >= TRACK_THRESHOLD ) {
2624 stripXY[1][i] = m_stripResults[0][i];
2625 stripEXY[1][i] = m_stripResults[1][i];
2626 }
2627 else {
2628 stripXY[1][i] = 0;
2629 stripEXY[1][i] = 0;
2630 }
2631 }
2632 m_geStripEff = new TGraphErrors(STRIP_MAX, stripXY[0], stripXY[1], stripEXY[0], stripEXY[1]);
2633 m_geStripEff->SetMarkerStyle(25);
2634 m_geStripEff->SetMarkerSize(0.5);
2635 m_cv[2] = new TCanvas("GoodStripEff","Strip efficiency", 100, 100, 800, 600);
2636 m_cv[2]->SetFillColor(0);
2637 m_cv[2]->SetBorderMode(0);
2638 m_geStripEff->Draw("AP");
2639 m_cv[2]->Write();
2640
2641 // Save histograms
2642 for(int i=0; i<B_LAY_NUM; i++ )
2643 {
2644 m_hHitMapBarrel_Lay[i]->Write();
2645
2646 if( i<E_LAY_NUM ) {
2647 m_hHitMapEndcap_Lay[0][i]->Write();
2648 m_hHitMapEndcap_Lay[1][i]->Write();
2649 }
2650 }
2651
2652 for(int i=0; i<B_SEG_NUM; i++)
2653 {
2654 m_hHitMapBarrel_Seg[i]->Write();
2655
2656 if( i< E_SEG_NUM ) {
2657 m_hHitMapEndcap_Seg[0][i]->Write();
2658 m_hHitMapEndcap_Seg[1][i]->Write();
2659 }
2660 }
2661 m_hTrackDistance->Fit("gaus");
2662 m_hTrackPosPhiDiff->Fit("gaus");
2663 m_hTrackPosThetaDiff->Fit("gaus");
2664 m_hTrackMomPhiDiff->Fit("gaus");
2665 m_hTrackMomThetaDiff->Fit("gaus");
2666
2667 m_hTrackDistance->Write();
2668 m_hTrackPosPhiDiff->Write();
2669 m_hTrackPosThetaDiff->Write();
2670 m_hTrackMomPhiDiff->Write();
2671 m_hTrackMomThetaDiff->Write();
2672
2673 for( int i=0; i<B_LAY_NUM; i++ ) m_hBarrelResDist[i]->Write();
2674 for( int i=0; i<E_LAY_NUM; i++ ) m_hEndcapResDist[i]->Write();
2675
2676 m_hBarrelResComp[0]->Write();
2677 m_hBarrelResComp[1]->Write();
2678 m_hEndcapResComp[0]->Write();
2679 m_hEndcapResComp[1]->Write();
2680
2681 if( m_usePad != 0 ) m_histArray->Write();
2682
2683 for( int i=0; i<BOX_MAX; i++ )
2684 {
2685 m_hStripFireMap[ i ]->Write();
2686 m_hStripExpHitMap[ i ]->Write();
2687 m_hStripEffHitMap[ i ]->Write();
2688 m_hStripNosHitMap[ i ]->Write();
2689 m_hStripEffMap[ i ]->Write();
2690 m_hStripNosRatioMap[ i ]->Write();
2691 }
2692 m_hStripFire->Write();
2693 m_hStripExpHit->Write();
2694 m_hStripEffHit->Write();
2695 m_hStripNosHit->Write();
2696 m_hStripEff->Write();
2697 m_hStripArea->Write();
2698 m_hStripNos->Write();
2699 m_hStripNosRatio->Write();
2700 log << MSG::INFO << "Save LV2 histograms done!" << endreq;
2701
2702 m_hBoxFire->Write();
2703 m_hBoxExpHit->Write();
2704 m_hBoxEffHit->Write();
2705 m_hBoxNosHit->Write();
2706 m_hBoxEff->Write();
2707 m_hBoxArea->Write();
2708 m_hBoxNos->Write();
2709 m_hBoxNosRatio->Write();
2710 log << MSG::INFO << "Save LV1 histograms done!" << endreq;
2711
2712 m_hBrLayerFire->Write();
2713 m_hEcLayerFire->Write();
2714 m_hLayerFire->Write();
2715 m_hLayerExpHit->Write();
2716 m_hLayerEffHit->Write();
2717 m_hLayerNosHit->Write();
2718 m_hLayerEff->Write();
2719 m_hLayerArea->Write();
2720 m_hLayerNos->Write();
2721 m_hLayerNosRatio->Write();
2722
2723 for( int i=0; i<LAYER_MAX; i++ ) m_hLayerCluster[i]->Write();
2724 for( int i=0; i<BOX_MAX; i++ ) m_hBoxCluster[i]->Write();
2725 m_hLayerClusterCmp->Write();
2726 m_hBoxClusterCmp->Write();
2727
2728 log << MSG::INFO << "Save histograms done!" << endreq;
2729
2730 // Save trees
2731 m_fLayerCoverage = 100*(double)m_fCalibLayerNum/LAYER_MAX;
2732 m_fBoxCoverage = 100*(double)m_fCalibBoxNum/BOX_MAX;
2733 m_fStripCoverage = 100*(double)m_fCalibStripNum/STRIP_MAX;
2734
2735 long digi_num, trk_num, eff_hit, nos_hit, exp_hit;
2736 m_tStatLog->Branch("digi_num", &digi_num, "digi_num/I");
2737 m_tStatLog->Branch("trk_num", &trk_num, "trk_num/I");
2738 m_tStatLog->Branch("eff_hit", &eff_hit, "eff_hit/I");
2739 m_tStatLog->Branch("nos_hit", &nos_hit, "nos_hit/I");
2740 m_tStatLog->Branch("exp_hit", &exp_hit, "exp_hit/I");
2741
2742 int stripMax;
2743 for( int i=0; i<PART_MAX; i++ )
2744 for( int j=0; j<((i==BRID)?B_SEG_NUM:E_SEG_NUM); j++ )
2745 for( int k=0; k<((i==BRID)?B_LAY_NUM:E_LAY_NUM); k++ )
2746 {
2747 stripMax = m_ptrIdTr->GetStripMax( i, j, k );
2748 for( int n=0; n<stripMax; n++ )
2749 {
2750 digi_num = m_record[i][j][k][n][0];
2751 trk_num = m_record[i][j][k][n][1];
2752 eff_hit = m_record[i][j][k][n][2];
2753 nos_hit = m_record[i][j][k][n][3];
2754 exp_hit = m_record[i][j][k][n][4];
2755 m_tStatLog->Fill();
2756 }
2757 }
2758
2759 m_jobFinish = clock();
2760 m_fTotalJobTime = (double)(m_jobFinish - m_jobStart)/CLOCKS_PER_SEC;
2761
2762 m_tJobLog->Fill();
2763 m_tJobLog->Write();
2764 m_tStatLog->Write();
2765 m_tLayConst->Write();
2766 m_tBoxConst->Write();
2767 m_tStrConst->Write();
2768
2769 // Close cluster output file
2770 if( m_fdata != NULL ) m_fdata->close();
2771
2772 return StatusCode::SUCCESS;
2773}
2774
2775StatusCode MucCalibMgr::EndRun()
2776{
2777 MsgStream log(msgSvc, "MucCalibMgr");
2778 log << MSG::INFO << endreq;
2779 log << MSG::INFO << "Total events : " << m_fTotalEvent << endreq;
2780 log << MSG::INFO << "Total digis : " << m_fTotalDigi << endreq;
2781 log << MSG::INFO << "Total rec hits : " << m_fTotalExpHit << endreq;
2782 log << MSG::INFO << "Total eff hits : " << m_fTotalEffHit << endreq;
2783 log << MSG::INFO << "Total inc hits : " << m_fTotalNosHit << endreq;
2784 log << MSG::INFO << "Total clusters : " << m_fTotalClstNum<< endreq;
2785
2786 log << MSG::INFO << "Strip calibrated percentage: "
2787 << 100*(double)m_fCalibStripNum/STRIP_MAX << "%" << endreq;
2788 log << MSG::INFO << "Box calibrated percentage: "
2789 << 100*(double)m_fCalibBoxNum/BOX_MAX << "%" << endreq;
2790 log << MSG::INFO << "Layer calibrated percentage: "
2791 << 100*(double)m_fCalibLayerNum/LAYER_MAX << "%" << endreq;
2792
2793 log << MSG::INFO << "Calibration run successfully" << endreq << endreq;
2794
2795 return StatusCode::SUCCESS;
2796}
2797
2798// END
Write
curve Write()
sprintf
sprintf(cut,"kal_costheta0_em>-0.93&&kal_costheta0_em<0.93&&kal_pxy0_em>=0.05+%d*0.1&&kal_pxy0_em<0.15+%d*0.1&&NGch>=2", j, j)
n
const Int_t n
Definition DataBase/tau_mode.c:65
phi2
Double_t phi2
Definition DataBase/tau_mode.c:7
phi1
Double_t phi1
Definition DataBase/tau_mode.c:7
e1
Double_t e1
Definition DataBase/tau_mode.c:7
e2
Double_t e2
Definition DataBase/tau_mode.c:7
PI
#define PI
Definition EvtD0ToKKpi0.cc:41
title
titledef title[20]
Definition EvtPhokharaDef.hh:61
STRIP_MAX
const int STRIP_MAX
Definition MucMappingAlg.h:24
mark_col
vector< MucMark * > mark_col
Definition MucMark.h:17
RecMucTrackCol
ObjectVector< RecMucTrack > RecMucTrackCol
Definition RecMucTrack.h:394
Segment
#define Segment
Definition TTrackBase.h:31
theta2
@ theta2
Definition TrkKalDeriv.h:24
theta1
@ theta1
Definition TrkKalDeriv.h:24
NULL
#define NULL
Definition VertexExtrapolate.h:14
MucCalibMgr::ClearHistoLV2
StatusCode ClearHistoLV2()
Definition MucCalibMgr.cxx:279
MucCalibMgr::InitResHisto
StatusCode InitResHisto()
Definition MucCalibMgr.cxx:869
MucCalibMgr::FillEffHit
StatusCode FillEffHit(int part, int segment, int layer, int strip)
Definition MucCalibMgr.cxx:2093
MucCalibMgr::DimuSelect
StatusCode DimuSelect()
Definition MucCalibMgr.cxx:1033
MucCalibMgr::InitOnlineHisto
StatusCode InitOnlineHisto()
Definition MucCalibMgr.cxx:602
MucCalibMgr::Init2DEffHisto
StatusCode Init2DEffHisto()
Definition MucCalibMgr.cxx:791
MucCalibMgr::SaveConst
StatusCode SaveConst()
Definition MucCalibMgr.cxx:2560
MucCalibMgr::FillExpHit
StatusCode FillExpHit(int part, int segment, int layer, int strip)
Definition MucCalibMgr.cxx:2080
MucCalibMgr::ClearResHisto
StatusCode ClearResHisto()
Definition MucCalibMgr.cxx:338
MucCalibMgr::ClearHistoLV0
StatusCode ClearHistoLV0()
Definition MucCalibMgr.cxx:245
MucCalibMgr::ClearClusterHisto
StatusCode ClearClusterHisto()
Definition MucCalibMgr.cxx:321
MucCalibMgr::AnalyseCluster
StatusCode AnalyseCluster()
Definition MucCalibMgr.cxx:2514
MucCalibMgr::ClearHistoLV1
StatusCode ClearHistoLV1()
Definition MucCalibMgr.cxx:263
MucCalibMgr::PadEff
StatusCode PadEff()
Definition MucCalibMgr.cxx:2478
MucCalibMgr::InitHistoLV1
StatusCode InitHistoLV1()
Definition MucCalibMgr.cxx:724
MucCalibMgr::InitHistoLV2
StatusCode InitHistoLV2()
Definition MucCalibMgr.cxx:740
MucCalibMgr::InitClusterHisto
StatusCode InitClusterHisto()
Definition MucCalibMgr.cxx:839
MucCalibMgr::FillCluster
StatusCode FillCluster(int part, int segment, int layer, int size)
Definition MucCalibMgr.cxx:2119
MucCalibMgr::ClearOnlineHisto
StatusCode ClearOnlineHisto()
Definition MucCalibMgr.cxx:224
MucCalibMgr::InitHisto
StatusCode InitHisto()
Definition MucCalibMgr.cxx:576
MucCalibMgr::FillEvent
StatusCode FillEvent()
Definition MucCalibMgr.cxx:1886
MucCalibMgr::EffAndNoiseLV2
StatusCode EffAndNoiseLV2()
Definition MucCalibMgr.cxx:2372
MucCalibMgr::ntupleSvc
INTupleSvc * ntupleSvc
Definition MucCalibMgr.h:100
MucCalibMgr::MucCalibMgr
MucCalibMgr(std::vector< double > jobInfo, std::vector< int > configInfo, std::string outputFile)
Definition MucCalibMgr.cxx:65
MucCalibMgr::InitHistoLV0
StatusCode InitHistoLV0()
Definition MucCalibMgr.cxx:705
MucCalibMgr::eventSvc
IDataProviderSvc * eventSvc
Definition MucCalibMgr.h:101
MucCalibMgr::EndRun
StatusCode EndRun()
Definition MucCalibMgr.cxx:2775
MucCalibMgr::ReadEvent
StatusCode ReadEvent()
Definition MucCalibMgr.cxx:1154
MucCalibMgr::EffAndNoiseLV1
StatusCode EffAndNoiseLV1()
Definition MucCalibMgr.cxx:2262
MucCalibMgr::AnalyseRes
StatusCode AnalyseRes()
Definition MucCalibMgr.cxx:2521
MucCalibMgr::msgSvc
IMessageSvc * msgSvc
Definition MucCalibMgr.h:99
MucCalibMgr::CheckEvent
StatusCode CheckEvent()
Definition MucCalibMgr.cxx:1831
MucCalibMgr::FillNosHit
StatusCode FillNosHit(int part, int segment, int layer, int strip)
Definition MucCalibMgr.cxx:2106
MucCalibMgr::EffAndNoiseLV0
StatusCode EffAndNoiseLV0()
Definition MucCalibMgr.cxx:2152
MucCalibMgr::InitConstTree
StatusCode InitConstTree()
Definition MucCalibMgr.cxx:897
MucCalibMgr::ClearConstTree
StatusCode ClearConstTree()
Definition MucCalibMgr.cxx:356
MucCalibMgr::InitArea
StatusCode InitArea()
Definition MucCalibMgr.cxx:992
MucCalibMgr::AnalyseEffAndNoise
StatusCode AnalyseEffAndNoise()
Definition MucCalibMgr.cxx:2134
MucCalibMgr::InitNtuple
StatusCode InitNtuple()
Definition MucCalibMgr.cxx:374
MucCalibMgr::Clear2DEffHisto
StatusCode Clear2DEffHisto()
Definition MucCalibMgr.cxx:305
MucCalibMgr::FillDigi
StatusCode FillDigi(int part, int segment, int layer, int strip)
Definition MucCalibMgr.cxx:2016
MucEntityCal::GetArea
double GetArea()
Definition MucEntityCal.cxx:95
MucID::barrel_ec
static int barrel_ec(const Identifier &id)
Values of different levels.
Definition MucID.cxx:41
MucID::layer
static int layer(const Identifier &id)
Definition MucID.cxx:61
MucID::channel
static int channel(const Identifier &id)
Definition MucID.cxx:71
MucID::segment
static int segment(const Identifier &id)
Definition MucID.cxx:51
MucIdTransform::GetStripId
int GetStripId(int part, int segment, int layer, int subid)
Definition Muc/MucCalibAlg/MucCalibAlg-00-03-00/src/MucIdTransform.cxx:67
MucIdTransform::SetBoxPos
bool SetBoxPos(int boxid, int *part, int *segment, int *layer)
Definition Muc/MucCalibAlg/MucCalibAlg-00-03-00/src/MucIdTransform.cxx:100
MucIdTransform::GetStripMax
int GetStripMax(int part, int segment, int layer)
Definition Muc/MucCalibAlg/MucCalibAlg-00-03-00/src/MucIdTransform.cxx:38
MucMark::IsInSegWith
bool IsInSegWith(MucMark &other)
Definition MucMark.cxx:144
MucMark::IsInCol
int IsInCol(int part, int segment, int layer, int strip, mark_col &aCol)
Definition MucMark.cxx:100
TofHitStatus::is_cluster
bool is_cluster() const
Definition TofHitStatus.h:25
TofHitStatus::setStatus
void setStatus(unsigned int status)
Definition TofHitStatus.cxx:134
Event
Definition Event.h:21