227 {
228
229 MsgStream log(
msgSvc(), name());
230 log << MSG::INFO << "in execute()" << endreq;
231
232 SmartDataPtr<Event::EventHeader> eventHeader(eventSvc(),"/Event/EventHeader");
233 m_run = eventHeader->runNumber();
234 m_event = eventHeader->eventNumber();
235 log << MSG::DEBUG <<"Run " << m_run << "\tEvent " << m_event <<endreq;
236
237 SmartDataPtr<DQAEvent::DQAEvent> dqaevt(eventSvc(), "/Event/DQATag");
238 if ( dqaevt ) {
239 log << MSG::INFO << "success get DQAEvent" << endreq;
240 } else {
241 log << MSG::ERROR << "Error accessing DQAEvent" << endreq;
242 return StatusCode::FAILURE;
243 }
244 log << MSG::DEBUG << "DQA event tag = " << dqaevt->EventTag() << endreq;
245
246 int part, segment, layer, strip;
247 char name[100];
248 bool isDimu = false;
249 if ( dqaevt->Dimu() ) isDimu = true;
250 log << MSG::INFO << "DQADimuTag:\t" << dqaevt->Dimu() << endreq;
251
252
253 log << MSG::INFO << "Retrieve digis" << endreq;
254 SmartDataPtr<MucDigiCol> mucDigiCol(eventSvc(),"/Event/Digi/MucDigiCol");
255 if(!mucDigiCol) {
256 log << MSG::FATAL << "Could not find MUC digi" << endreq;
257 return( StatusCode::FAILURE);
258 }
259
261 MucDigiCol::iterator digiIter = mucDigiCol->begin();
262 int eventDigi = 0;
263 for ( int digiId =0; digiIter != mucDigiCol->end(); digiIter++, digiId++ )
264 {
265 mucId = (*digiIter)->identify();
268 eventDigi ++;
269
270
272 m_digiCol.push_back( aMark );
273 m_segDigiCol[part][segment].push_back( aMark );
274 m_recordAll[part][segment][layer][strip][0] ++;
275 if( isDimu ) m_recordDimu[part][segment][layer][strip][0] ++;
276 }
277 log << MSG::INFO << "Total digis in this event: " << eventDigi << endreq;
278 if( eventDigi > 500) cout <<"Run:\t"<< m_run <<"\tEvent:\t"<< m_event << "\tdigits inflation:\t" << eventDigi << endreq;
279
280
281 m_clusterCol = (*m_ptrMucMark).CreateClusterCol(1, m_digiCol );
282 log << MSG::INFO << "Total clusters in this event: " << m_clusterCol.size() << endreq;
283
284
285
287 SmartDataPtr<RecMucTrackCol> mucTrackCol(eventSvc(),"/Event/Recon/RecMucTrackCol");
288 if (!mucTrackCol) {
289 log << MSG::FATAL << "Could not find RecMucTrackCol" << endreq;
290 return( StatusCode::FAILURE);
291 }
292
294 if (aRecMucTrackCol->size() < 1) {
295 log << MSG::INFO << "No MUC tracks in this event" << endreq;
296 return StatusCode::SUCCESS;
297 }
298 log << MSG::INFO << "Total tracks of this event: " << aRecMucTrackCol->size() << endreq;
299
300 int trackHitNum, expectedHitNum, segNum, lastLayerBR, lastLayerEC;
301 int layerPassNum[3], passMax[TRACK_SEG_MAX][2];
302 bool firedLay[TRACK_SEG_MAX][LAYER_MAX];
303 double costheta, phi;
304 TH1* htmp(0);
305 trackHitNum = expectedHitNum = segNum = lastLayerBR = lastLayerEC = 0;
306 layerPassNum[0] = layerPassNum[1] = layerPassNum[2] = 0;
307 for( int segi=0; segi<TRACK_SEG_MAX; segi++ ) {
308 passMax[segi][0] = passMax[segi][1] = 0;
309 for( int layi=0; layi<LAYER_MAX; layi++ ) firedLay[segi][layi] = 0;
310 }
311
312 vector<MucRecHit*> mucRawHitCol;
313 vector<MucRecHit*> mucExpHitCol;
314 RecMucTrackCol::iterator trackIter = mucTrackCol->begin();
315 for (int trackId = 0; trackIter != mucTrackCol->end(); trackIter++, trackId++)
316 {
317 trackHitNum = (*trackIter)->numHits();
318
319 if( trackHitNum == 0 ) {
320 log << MSG::INFO << "Track " << trackId << " no hits" << endreq;
321 continue;
322 }
323
324 lastLayerBR = (*trackIter)->brLastLayer();
325 lastLayerEC = (*trackIter)->ecLastLayer();
326
327 CLHEP::Hep3Vector a3Vector((*trackIter)->xPos(),(*trackIter)->yPos(),(*trackIter)->zPos());
328 costheta = a3Vector.cosTheta();
329 phi = a3Vector.phi();
330 if(m_thsvc->getHist("/DQAHist/MUC/Costheta_All", htmp).isSuccess()) {
331 htmp->Fill( costheta );
332 } else {
333 log << MSG::ERROR << "Fail to retrieve Costheta_All" << endreq;
334 }
335 if(m_thsvc->getHist("/DQAHist/MUC/Phi_All", htmp).isSuccess()) {
336 htmp->Fill( phi );
337 } else {
338 log << MSG::ERROR << "Fail to retrieve Phi_All" << endreq;
339 }
340
341 if( isDimu )
342 {
343 if(m_thsvc->getHist("/DQAHist/MUC/Costheta_Dimu", htmp).isSuccess()) {
344 htmp->Fill( costheta );
345 } else {
346 log << MSG::ERROR << "Fail to retrieve Costheta_Dimu" << endreq;
347 }
348 if(m_thsvc->getHist("/DQAHist/MUC/Phi_Dimu", htmp).isSuccess()) {
349 htmp->Fill( phi );
350 } else {
351 log << MSG::ERROR << "Fail to retrieve Phi_Dimu" << endreq;
352 }
353 }
354 log << MSG::INFO << "Fill costheta and phi:\t" << costheta << "\t" << phi << endreq;
355
358
359
360
361 mucExpHitCol = (*trackIter)->GetExpectedHits();
362
363 expectedHitNum += mucExpHitCol.size();
364 for(unsigned int hitId = 0; hitId < mucExpHitCol.size(); hitId++)
365 {
366 pMucRawHit = mucExpHitCol[ hitId ];
367 part = pMucRawHit->
Part(); segment = pMucRawHit->
Seg();
368 layer = pMucRawHit->
Gap(); strip = pMucRawHit->
Strip();
369
370 MucMark* currentMark =
new MucMark( part, segment, layer, strip );
371 m_expHitCol.push_back( currentMark );
372
373
374 int isInPos = -1;
375 bool isInEffWindow = false;
376 isInPos = currentMark->
IsInCol( m_segDigiCol[part][segment] );
377
378
379 if( part == BRID && (layer-lastLayerBR>1) ) continue;
380 if( part != BRID && (layer-lastLayerEC>1) ) continue;
381
382
383 if( part==BRID && layer==0 && (strip<2 || strip>45) )
384 {
385 if( isInPos != -1)
386 {
387 m_recordAll[part][segment][layer][strip][2] ++;
388 m_recordAll[part][segment][layer][strip][1] ++;
389 m_effHitCol.push_back( m_segDigiCol[part][segment][isInPos] );
390
391 if(isDimu) {
392 m_recordDimu[part][segment][layer][strip][2] ++;
393 m_recordDimu[part][segment][layer][strip][1] ++;
394 }
395 }
396 else {
397 m_recordAll[part][segment][layer][strip][1] ++;
398 if(isDimu) m_recordDimu[part][segment][layer][strip][1] ++;
399 }
400 continue;
401 }
402
403 if( isInPos != -1 )
404 {
405 m_recordAll[part][segment][layer][strip][2] ++;
406 m_recordAll[part][segment][layer][strip][1] ++;
407 m_effHitCol.push_back( m_segDigiCol[part][segment][isInPos] );
408
409 if(isDimu) {
410 m_recordDimu[part][segment][layer][strip][2] ++;
411 m_recordDimu[part][segment][layer][strip][1] ++;
412 }
413
414 continue;
415 }
416 else for(int tempStrip=0, hiti=-m_effWindow; hiti<=m_effWindow; hiti++ )
417 {
418 if( hiti == 0 ) continue;
419 tempStrip = strip + hiti;
420 if( tempStrip < 0 || tempStrip > m_ptrIdTr->
GetStripMax(part,segment,layer) )
continue;
421
422 isInPos = m_ptrMucMark->
IsInCol( part, segment, layer, tempStrip, m_segDigiCol[part][segment] );
423 if( isInPos != -1 )
424 {
425 m_recordAll[part][segment][layer][tempStrip][2] ++;
426 m_recordAll[part][segment][layer][tempStrip][1] ++;
427 m_effHitCol.push_back( m_segDigiCol[part][segment][isInPos] );
428
429 if(isDimu) {
430 m_recordDimu[part][segment][layer][tempStrip][2] ++;
431 m_recordDimu[part][segment][layer][tempStrip][1] ++;
432 }
433
434 isInEffWindow = true;
435 }
436 }
437
438 if( isInEffWindow ) { continue; }
439 else {
440 m_recordAll[part][segment][layer][strip][1] ++;
441 if(isDimu) m_recordDimu[part][segment][layer][strip][1] ++;
442 }
443
444 }
445
446
447 log << MSG::INFO << "Fill residual" << endreq;
448 vector<float> m_lineResCol = (*trackIter)->getDistHits();
449
450 mucRawHitCol = (*trackIter)->GetHits();
451
452
453 if( trackHitNum > 4 && m_lineResCol[0] != -99)
454 {
455
456 bool firedFlag[
PART_MAX][LAYER_MAX][2];
457 for(
int iprt=0; iprt<
PART_MAX; iprt++)
458 for(int jlay=0; jlay<LAYER_MAX; jlay++)
459 firedFlag[iprt][jlay][0] = firedFlag[iprt][jlay][1] = false;
460
461 for(unsigned int hitId = 0; hitId < mucExpHitCol.size(); hitId++)
462 {
463 pMucExpHit = mucExpHitCol[ hitId ];
464 part = pMucExpHit->
Part(); segment = pMucExpHit->
Seg(); layer = pMucExpHit->
Gap();
465 firedFlag[part][layer][0] = true;
466 }
467
468 for(unsigned int hitId = 0; hitId < mucRawHitCol.size(); hitId++)
469 {
470 pMucRawHit = mucRawHitCol[ hitId ];
471 part = pMucRawHit->
Part(); segment = pMucRawHit->
Seg(); layer = pMucRawHit->
Gap();
472
473
474 if( firedFlag[part][layer][0] == true && firedFlag[part][layer][1] == false )
475 {
476 if( part == BRID )
477 {
478 sprintf( name, "/DQAHist/MUC/BrResDist_All_L%d", layer );
479 if(m_thsvc->getHist(name, htmp).isSuccess()) {
480 htmp->Fill( m_lineResCol[hitId] );
481 } else {
482 log << MSG::ERROR << "Fail to retrieve " << name << endreq;
483 }
484
485
486 if( isDimu )
487 {
488 sprintf( name, "/DQAHist/MUC/BrResDist_Dimu_L%d", layer );
489 if(m_thsvc->getHist(name, htmp).isSuccess()) {
490 htmp->Fill( m_lineResCol[hitId] );
491 } else {
492 log << MSG::ERROR << "Fail to retrieve " << name << endreq;
493 }
494
495 }
496 } else {
497 sprintf( name, "/DQAHist/MUC/EcResDist_All_L%d", layer );
498 if(m_thsvc->getHist(name, htmp).isSuccess()) {
499 htmp->Fill( m_lineResCol[hitId] );
500 } else {
501 log << MSG::ERROR << "Fail to retrieve " << name << endreq;
502 }
503
504
505 if( isDimu )
506 {
507 sprintf( name, "/DQAHist/MUC/EcResDist_Dimu_L%d", layer );
508 if(m_thsvc->getHist(name, htmp).isSuccess()) {
509 htmp->Fill( m_lineResCol[hitId] );
510 } else {
511 log << MSG::ERROR << "Fail to retrieve " << name << endreq;
512 }
513
514 }
515 }
516 }
517
518 firedFlag[part][layer][1] = true;
519 }
520 }
521
522 mucRawHitCol.clear();
523 mucExpHitCol.clear();
524
525 }
526
527
528
529 log << MSG::INFO << "Searching inc/noise hits" << endreq;
530 bool isNosHit;
531 bool hasEffHit;
532 for( unsigned int i=0; i < m_digiCol.size(); i++ )
533 {
534 isNosHit = true;
535
536 if( m_digiCol[i]->IsInCol( m_effHitCol )!=-1 ) continue;
537 else
538 {
539 for( unsigned int j=0; j < m_clusterCol.size(); j++ )
540 {
541 hasEffHit = false;
542 for( unsigned int k=0; k<m_clusterCol[j].size(); k++)
543 {
544 if( m_clusterCol[j][k]->IsInCol(m_effHitCol) != -1)
545 {
546 hasEffHit = true;
547 break;
548 }
549 }
550
551 if ( hasEffHit && (m_digiCol[i]->IsInCol( m_clusterCol[j] ) != -1) ) {
552 isNosHit = false;
553 break;
554 }
555 }
556
557 if( isNosHit ) {
558 part = (*m_digiCol[i]).Part(); segment = (*m_digiCol[i]).
Segment();
559 layer = (*m_digiCol[i]).Layer(); strip = (*m_digiCol[i]).Strip();
560
561
562 m_recordAll[part][segment][layer][strip][3] ++;
563 if( isDimu ) m_recordDimu[part][segment][layer][strip][3] ++;
564 m_nosHitCol.push_back(m_digiCol[i]);
565 }
566 }
567 }
568
569
570 DQA_MUC::FillHistograms( isDimu );
571
572
573 for( unsigned int i=0; i<m_digiCol.size(); i++ ) {
574 if( m_digiCol[i] != NULL ) delete m_digiCol[i];
575 }
576
577 for( unsigned int i=0; i<m_expHitCol.size(); i++ ) {
578 if( m_expHitCol[i] != NULL ) delete m_expHitCol[i];
579 }
580
581 if( m_effHitCol.size() != 0 ) m_effHitCol.clear();
582 if( m_expHitCol.size() != 0 ) m_expHitCol.clear();
583 if( m_nosHitCol.size() != 0 ) m_nosHitCol.clear();
584 if( m_digiCol.size() != 0 ) m_digiCol.clear();
585 if( m_clusterCol.size() != 0 ) m_clusterCol.clear();
586
588 for( int j=0; j<SEGMENT_MAX; j++ ) {
589 if( m_segDigiCol[i][j].size() != 0 ) m_segDigiCol[i][j].clear();
590 }
591 }
592
593 return StatusCode::SUCCESS;
594}
ObjectVector< RecMucTrack > RecMucTrackCol
static int barrel_ec(const Identifier &id)
Values of different levels.
static int layer(const Identifier &id)
static int channel(const Identifier &id)
static int segment(const Identifier &id)
int IsInCol(int part, int segment, int layer, int strip, mark_col &aCol)
int Part() const
Get Part.
int Strip() const
Get Strip.
_EXTERN_ std::string EvtRecEvent