Changeset 550 in svn
- Timestamp:
- Feb 22, 2010, 7:38:15 PM (15 years ago)
- Location:
- trunk
- Files:
-
- 6 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/Delphes.cpp
r547 r550 628 628 elementElec->IsolFlag = DET->Isolation(electron[i],TrackCentral,DET->ISOL_PT,DET->ISOL_Cone,ptisoEl); 629 629 elementElec->IsolPt = ptisoEl; 630 DET->BinEtaPhi(elementElec->PhiCalo,elementElec->EtaCalo,iPhiEl,iEtaEl); 630 631 cout << "\nC'est un electron\n"; 632 int electron_tower_index = DET->BinEtaPhi(elementElec->PhiCalo,elementElec->EtaCalo,iPhiEl,iEtaEl); 631 633 D_CaloTower calElec(list_of_active_towers.getElement(iEtaEl,iPhiEl)); 632 634 elementElec->EHoverEE = calElec.getEhad()/calElec.getEem(); 635 elementElec->EtRatio = DET->CaloIsolation(muon[i], list_of_active_towers,iPhiEl,iEtaEl,electron_tower_index); 633 636 } 634 637 … … 644 647 elementMu->IsolFlag = DET->Isolation(muon[i],TrackCentral,DET->ISOL_PT,DET->ISOL_Cone,ptisoMu); 645 648 elementMu->IsolPt = ptisoMu; 646 DET->BinEtaPhi(elementMu->PhiCalo,elementMu->EtaCalo,iPhiMu,iEtaMu); 649 650 cout << "\nC'est un muon\n"; 651 // binning of the muon position in calo 652 int muon_tower_index = DET->BinEtaPhi(elementMu->PhiCalo,elementMu->EtaCalo,iPhiMu,iEtaMu); 647 653 D_CaloTower calMuon(list_of_active_towers.getElement(iEtaMu,iPhiMu)); 648 654 if( calMuon.getEem() !=0 ) elementMu->EHoverEE = calMuon.getEhad()/calMuon.getEem(); 649 655 else elementMu->EHoverEE = UNDEFINED; 650 elementMu->EtRatio = DET->CaloIsolation(muon[i], list_of_active_towers,iPhiMu,iEtaMu );656 elementMu->EtRatio = DET->CaloIsolation(muon[i], list_of_active_towers,iPhiMu,iEtaMu,muon_tower_index); 651 657 } 652 658 -
trunk/Utilities/ExRootAnalysis/interface/BlockClasses.h
r547 r550 262 262 float PhiCalo; // particle azimuthal angle in rad when entering the calo 263 263 float EHoverEE; 264 float EtRatio; 264 265 265 266 void SetEtaPhiCalo(const float eta, const float phi) {EtaCalo=eta; PhiCalo=phi;}; -
trunk/Utilities/ExRootAnalysis/src/BlockClasses.cc
r539 r550 29 29 ** ** 30 30 ***********************************************************************/ 31 31 32 32 33 #include "BlockClasses.h" -
trunk/interface/SmearUtil.h
r531 r550 278 278 279 279 /// Lepton isolation based on calorimetry (optional. Default: off) 280 float CaloIsolation(const D_Particle& part, const D_CaloTowerList & towers, const float iPhi, const float iEta );280 float CaloIsolation(const D_Particle& part, const D_CaloTowerList & towers, const float iPhi, const float iEta, const int iTower); 281 281 282 282 //********************* returns a segmented value for eta and phi, for calo towers ***** 283 voidBinEtaPhi(const float phi, const float eta, float& iPhi, float& iEta);283 int BinEtaPhi(const float phi, const float eta, float& iPhi, float& iEta); 284 284 285 285 void setNames(const string& list, const string& det, const string& trig); -
trunk/src/CaloUtil.cc
r443 r550 237 237 return _calotowers[i]; 238 238 } 239 cout << "Xav!!!!! not found in getElement" << endl; 239 240 return dummyTower; // if not found 240 241 } … … 247 248 return _calotowers[i]; 248 249 } 250 cout << "Xav!!!!! not found in getElement const" << endl; 249 251 return dummyTower; // if not found 250 252 } -
trunk/src/SmearUtil.cc
r537 r550 1440 1440 1441 1441 // ******* Calorimetric isolation 1442 float RESOLution::CaloIsolation(const D_Particle& part, const D_CaloTowerList & towers, const float iPhi, const float iEta) { 1443 // etrat, which is a percentage between 00 and 99. It is the ratio of the transverse energy 1444 // in a 3Ã3 grid surrounding the muon to the pT of the muon. For well-isolated muons, both ptiso and etrat will be small. 1442 float RESOLution::CaloIsolation(const D_Particle& part, const D_CaloTowerList & towers, const float iPhi, const float iEta, const int iTower) { 1443 // iPhi and iEta are the coordinates of the center of the tower crossed by the particle 1444 // iPhi is in [-pi ; pi] and iEta is in [-etamax ; etamax] 1445 // iTower is the index of the tower, in [0, n_tower]. iTower points only towers in positive range 1446 1447 // etrat, which is a percentage between 00 and 99. It is the ratio of the transverse energy 1448 // in a NÃN grid surrounding the muon to the pT of the muon. For well-isolated muons, both ptiso and etrat will be small. 1449 // available parameters: ISOL_Calo_ET , ISOL_Calo_Grid 1450 1445 1451 if(ISOL_Calo_ET>1E10) return UNDEFINED; // avoid doing anything unreasonable... 1446 1452 float et_sum=0; 1447 // available parameters: ISOL_Calo_ET , ISOL_Calo_Grid1448 // Get the EtaCalo/PhiCalo of the muon ;1449 // transform it into iEta/iPhi to get the towers, and their neighbourh (i-1, i-2, etc)1450 1453 1451 1454 unsigned int N = ISOL_Calo_Grid; 1452 int index= iUNDEFINED; // index of the central tower of the grid in TOWER_eta_edges[.]; 1453 // !! TOWER_eta_edges is only with eta>0 1454 // finds the index of the central tower of the NxN grid 1455 for (unsigned int i=1; i< TOWER_number+1; i++) { 1456 if(fabs(iEta) >= TOWER_eta_edges[i-1] && fabs(iEta) < TOWER_eta_edges[i]) { 1457 index = i-1; 1458 break; 1459 } 1455 if(N%2 ==0 || N==0) { cout << "Warning: ISOL_Calo_Grid must be odd. ISOL_Calo_Grid = 3 will be used\n"; N=3;} 1456 int index= iTower; // index of the central tower of the grid in TOWER_eta_edges[.]; 1457 // !! TOWER_eta_edges is only with eta>0 1458 1459 cout << "iEta= " << iEta << "\tiPhi= " << iPhi << endl; 1460 //cout << "la particule est dans la " << index << "eme cellule, entre " << TOWER_eta_edges[index] << " et " << TOWER_eta_edges[index+1] << endl; 1461 //float centre = (iEta>0) ? ((TOWER_eta_edges[index]+TOWER_eta_edges[index+1])/2.0) : -((TOWER_eta_edges[index]+TOWER_eta_edges[index+1])/2.0); 1462 //cout << "le centre de cette tour est " << centre << endl; 1463 1464 if(index != iUNDEFINED) { 1465 int range = (int) (N-1)/2; // the (int) conversion is needed, as N is "unsigned int" 1466 for(int i= -range; i<= range; i++) { // shift of index in eta: e.g.: i = -1, 0, 1 if N=3; 1467 // goal : to identify the center of each cell in the NxN grid. 1468 // the eta/phi coordinates of each center will be in (eta_ith_tower,phi_ith_tower) 1469 1470 unsigned int i_eta = (index+i>=0) ? index + i : (int) -1*(index + i +1); // i_eta is the index in [0; index_eta_max] 1471 //cout << i_eta << "\t"; 1472 float dphi = TOWER_dphi[i_eta]*pi/180.; // in rad // size in phi of the cells for this eta 1473 float eta_ith_tower = (TOWER_eta_edges[i_eta]+TOWER_eta_edges[i_eta+1])/2.0; 1474 if(iEta<0) eta_ith_tower *= -1; // needed if the central tower is in negative region 1475 if(index+i<0) eta_ith_tower *= -1; // needed if the central part is crossed by the grid 1476 //cout << "pour eta_ith_tower=" << eta_ith_tower << ", on va dphi = " << dphi << endl; 1477 // !!! at this point, eta_ith_tower is the value in eta of the center of the current calo cell 1478 1479 for(int j= -range; j<= range; j++) { // shift of index in phi. e.g.: j = -1, 0, 1 if N=3; 1480 float phi_ith_tower = iPhi + j*dphi; // in radian 1481 //cout << "eta_ith_tower= " << eta_ith_tower << "\tphi_ith_tower= " << phi_ith_tower << "\tj= " << j << "\tdphi=" << dphi << endl; 1482 if(phi_ith_tower > pi) phi_ith_tower -= 2.*pi; 1483 else if (phi_ith_tower < -pi) phi_ith_tower += 2.*pi; 1484 1485 D_CaloTower calMuon(towers.getElement(eta_ith_tower,phi_ith_tower)); 1486 cout << "eta_ith_tower= " << eta_ith_tower << "\tphi_ith_tower= " << phi_ith_tower <<"\t" 1487 << "calMuon.getEta= " << calMuon.getEta() << "\tcalMuon.getPhi()= " << calMuon.getPhi() <<"\t"; 1488 1489 if(calMuon.getEta() != UNDEFINED ) 1490 if(calMuon.getET() > ISOL_Calo_ET) { 1491 et_sum += calMuon.getET(); 1492 //cout << "eta/phi/Et = " << eta_ith_tower << "\t" << phi_ith_tower << "\t" << calMuon.getET() << " GeV" << endl; 1493 cout << calMuon.getET() << " GeV"; 1494 } 1495 //else cout << "eta/phi = " << eta_ith_tower << "\t" << phi_ith_tower << "\tnot active\n"; 1496 else cout << phi_ith_tower << "\tnot active"; 1497 cout << endl; 1498 } 1499 } 1500 } // undefined index 1501 else { 1502 cout << "out of range" << endl; 1503 if (CEN_max_mu < CEN_max_calo_fwd) 1504 cout << "** ERROR in RESOLution::CaloIsolation: 'muon'-tower not found! **" << endl; 1460 1505 } 1461 if(index != iUNDEFINED) {1462 // finds the size in phi of the cells for this eta1463 float dphi = TOWER_dphi[index]*pi/180.; // in rad1464 1465 //cout << "Grid " << " ----------\n";1466 for (unsigned int i_eta=0; i_eta<N; i_eta++) {1467 unsigned int real_index = (iEta>0) ? index+i_eta-(N-1)/2 : index+1+i_eta-(N-1)/2 ;1468 float eta_ith_tower = TOWER_eta_edges[real_index];1469 if(iEta<0) eta_ith_tower *= -1;1470 1471 for (unsigned int i_phi=0; i_phi<N; i_phi++) {1472 float phi_ith_tower = iPhi + (float)(i_phi - (N-1.)/2.)*dphi;1473 D_CaloTower calMuon(towers.getElement(eta_ith_tower,phi_ith_tower));1474 if(calMuon.getEta() != UNDEFINED && calMuon.getE() > ISOL_Calo_ET) {1475 et_sum += calMuon.getE();1476 //cout << "eta/phi = " << eta_ith_tower << "\t" << phi_ith_tower << "\t" << calMuon.getE() << " GeV" << endl;1477 }1478 //else cout << "eta/phi = " << eta_ith_tower << "\t" << phi_ith_tower << "\tnot active\n";1479 }1480 1481 } // NxN grid1482 //cout << "-----------" << etrat << endl;1483 }1484 else if (CEN_max_mu < CEN_max_calo_fwd)1485 cout << "** ERROR in RESOLution::CaloIsolation: 'muon'-tower not found! **" << endl;1486 1506 // should never happen ! this would be a bug 1487 1507 //cout << "etrat = " << et_sum << "\t Pt=" << part.Pt() << endl; … … 1490 1510 // which is Pt(muon) / sum(ET) 1491 1511 float etrat = 0.; 1492 if(et_sum==0) etrat = 99.;1493 else if(et_sum>0) etrat = 100*part.Pt()/et_sum;1512 /* if(et_sum==0) etrat = 99.; 1513 //else if(et_sum>0) etrat = 100*part.Pt()/et_sum; 1494 1514 if(etrat<0) cout << "Error: negative etrat in CaloIsolation (" << etrat <<")\n"; 1495 1515 //else if(etrat>99) cout << "Error: etrat should be in [0;99] in CaloIsolation (" << etrat <<")\n"; 1496 1516 if(etrat>99) etrat = 99; 1517 */ 1518 etrat = et_sum/part.Pt(); 1497 1519 return etrat; 1498 1520 } … … 1500 1522 1501 1523 //********** returns a segmented value for eta and phi, for calo towers ***** 1502 voidRESOLution::BinEtaPhi(const float phi, const float eta, float& iPhi, float& iEta){1524 int RESOLution::BinEtaPhi(const float phi, const float eta, float& iPhi, float& iEta){ 1503 1525 iEta = UNDEFINED; 1504 1526 int index= iUNDEFINED; … … 1510 1532 } 1511 1533 } 1512 if(index==UNDEFINED) return;1513 1534 iPhi = UNDEFINED; 1535 if(index==iUNDEFINED) return index; 1536 1514 1537 float dphi = TOWER_dphi[index]*pi/180.; 1515 1538 for (unsigned int i=1; i < 360/TOWER_dphi[index]; i++ ) { … … 1522 1545 } 1523 1546 if (phi > pi-dphi) iPhi = pi-dphi; 1547 return index; 1524 1548 } 1525 1549
Note:
See TracChangeset
for help on using the changeset viewer.