ECAL TB 2006 calibration with matrix inversion technique. More...

#include <InvRingCalib.h>

Inheritance diagram for InvRingCalib:

Public Member Functions
void	beginOfJob ()
	BeginOfJob.
Status	duringLoop (const edm::Event &, const edm::EventSetup &)
	duringLoop
void	endOfJob ()
Status	endOfLoop (const edm::EventSetup &, unsigned int iCounter)
	InvRingCalib (const edm::ParameterSet &)
	ctor
void	startingNewLoop (unsigned int)
	startingNewLoop
	~InvRingCalib ()
	dtor
Private Member Functions
int	EBRegId (const int)
	Gives back in which region you are:
void	EBRegionDef ()
	Defines the regions in the barrel.
int	EBRegionNum () const
	Number of regions in EB.
int	EERegId (int)
	gives back in which region of the endcap you are.
void	EERegionDef ()
	Defines the regions in the endcap.
int	EERegionNum () const
	The number of regions in EE.
void	EERingDef (const edm::EventSetup &)
	Defines the rins in the endcap.
void	RegPrepare ()
	Prepares the EB regions;.
Private Attributes
bool	isfirstcall_
edm::InputTag	m_barrelAlCa
	EcalBarrel Input Collection name.
std::vector< DetId >	m_barrelCells
	geometry things used all over the file
EcalIntercalibConstantMap	m_barrelMap
std::map< int, int >	m_cellPhi
std::map< int, GlobalPoint >	m_cellPos
	position of the cell, borders, coords etc...
std::string	m_EBcoeffFile
	coeffs filenames
std::string	m_EEcoeffFile
int	m_EEZone
	endcap zone to be calibrated
edm::InputTag	m_ElectronLabel
	To take the electrons.
edm::InputTag	m_endcapAlCa
	EcalEndcap Input Collection name.
std::vector< DetId >	m_endcapCells
EcalIntercalibConstantMap	m_endcapMap
int	m_endRing
int	m_etaEnd
	eta end of the zone of interest
int	m_etaStart
	eta start of the zone of interest
int	m_etaWidth
	eta size of the regions
std::vector< VEcalCalibBlock * >	m_IMACalibBlocks
	single blocks calibrators
unsigned int	m_loops
	LP sets the number of loops to do.
VFillMap *	m_MapFiller
	The class that fills the map!
std::string	m_mapFillerType
double	m_maxCoeff
	maximum coefficient accepted (RAW)
double	m_maxEnergyPerCrystal
	maximum energy per crystal cut
int	m_maxSelectedNumPerRing
	maximum number of events per Ring
double	m_minCoeff
	minimum coefficient accepted (RAW)
double	m_minEnergyPerCrystal
	minimum energy per crystal cut
int	m_recoWindowSidex
	reconstruction window size
int	m_recoWindowSidey
std::map< int, int >	m_Reg
	EB regions vs. eta index.
std::map< int, int >	m_RinginRegion
	association map between raw detIds and the number of the ring inside the region
std::map< int, int >	m_RingNumOfHits
	number of events already read per Ring
int	m_startRing
	LP define the EE region to calibrate.
int	m_usingBlockSolver
	to exclude the blocksolver
std::map< int, int >	m_xtalRegionId
	association map between raw detIds and Region
std::map< int, int >	m_xtalRing
	association map between Raw detIds and Rings

Detailed Description

ECAL TB 2006 calibration with matrix inversion technique.

Date:: 2010/01/18 21:31:47

Revision:: 1.5

Author:

Definition at line 29 of file InvRingCalib.h.

Constructor & Destructor Documentation

InvRingCalib::InvRingCalib ( const edm::ParameterSet & iConfig ) [explicit]

ctor

Definition at line 30 of file InvRingCalib.cc.

References abs, algorithm(), EBRegionNum(), EERegionNum(), cmsRelvalreport::exit, edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), i, m_EEZone, m_endRing, m_etaEnd, m_etaStart, m_etaWidth, m_IMACalibBlocks, m_loops, m_mapFillerType, m_startRing, and AlCaHLTBitMon_QueryRunRegistry::string.

                                                          :
  m_barrelAlCa (iConfig.getParameter<edm::InputTag>("barrelAlca")),
  m_endcapAlCa (iConfig.getParameter<edm::InputTag>("endcapAlca")) ,
  m_ElectronLabel (iConfig.getParameter<edm::InputTag>("ElectronLabel")),
  m_recoWindowSidex (iConfig.getParameter<int>("recoWindowSidex")),
  m_recoWindowSidey (iConfig.getParameter<int>("recoWindowSidey")),
  m_minEnergyPerCrystal (iConfig.getParameter<double>("minEnergyPerCrystal")) ,
  m_maxEnergyPerCrystal (iConfig.getParameter<double>("maxEnergyPerCrystal")) ,
  m_etaStart (iConfig.getParameter<int>("etaStart")),
  m_etaEnd (iConfig.getParameter<int>("etaEnd")),
  m_etaWidth (iConfig.getParameter<int>("etaWidth")),
  m_maxSelectedNumPerRing (iConfig.getParameter<int>("maxNumPerRing")),
  m_minCoeff (iConfig.getParameter<double>("minCoeff")),
  m_maxCoeff (iConfig.getParameter<double>("maxCoeff")),
  m_usingBlockSolver(iConfig.getParameter<int>("usingBlockSolver")),
  m_startRing (iConfig.getParameter<int>("startRing")),
  m_endRing (iConfig.getParameter<int>("endRing")),
  m_EBcoeffFile (iConfig.getParameter<std::string>("EBcoeffs")),
  m_EEcoeffFile (iConfig.getParameter<std::string>("EEcoeffs")),
  m_EEZone (iConfig.getParameter<int>("EEZone"))
{
  //controls if the parameters inputed are correct
  if ((m_etaEnd*m_etaStart)>0)
    assert (!((m_etaEnd - m_etaStart )%m_etaWidth)); 
  if ((m_etaEnd*m_etaStart)<0)
    assert (!((m_etaEnd - m_etaStart-1 )%m_etaWidth)); 

  assert (m_etaStart >=-85 && m_etaStart <= 86);
  assert (m_etaEnd >= m_etaStart && m_etaEnd <= 86);
  assert (m_startRing>-1 && m_startRing<= 40);
  assert (m_endRing>=m_startRing && m_endRing<=40);

  assert (!((m_endRing - m_startRing)%m_etaWidth));
  assert (( abs(m_EEZone)<=1));
  
  m_loops = (unsigned int) iConfig.getParameter<int>("loops")- 1;
  //LP CalibBlock vector instantiation
  edm::LogInfo ("IML") << "[InvRingCalib][ctor] Calib Block" ;
  std::string algorithm = iConfig.getParameter<std::string> ("algorithm") ;
  m_mapFillerType = iConfig.getParameter<std::string> ("FillType");
  int eventWeight = iConfig.getUntrackedParameter<int> ("L3EventWeight",1) ;
  
  for (int i = 0 ; i < EBRegionNum () ; ++i)
   {
  if (algorithm == "IMA")
        m_IMACalibBlocks.push_back (
            new IMACalibBlock (m_etaWidth)
          ) ; 
      else if (algorithm == "L3")
        m_IMACalibBlocks.push_back (
            new L3CalibBlock (m_etaWidth, eventWeight)
          ) ; 
      else
        {
          edm::LogError ("building") << algorithm 
                          << " is not a valid calibration algorithm" ;
          exit (1) ;    
        }      
   }   
  int EEBlocks = 0 ;
  if (m_EEZone == 0) EEBlocks = 2 * EERegionNum () ;
  if (m_EEZone == 1 || m_EEZone == -1) EEBlocks = EERegionNum () ;

  for (int i = 0; i < EEBlocks ; ++i)
   {
  
  if (algorithm == "IMA")
        m_IMACalibBlocks.push_back (
            new IMACalibBlock (m_etaWidth)
          ) ; 
      else if (algorithm == "L3")
        m_IMACalibBlocks.push_back (
            new L3CalibBlock (m_etaWidth, eventWeight)
          ) ; 
      else
        {
          edm::LogError ("building") << algorithm 
                          << " is not a valid calibration algorithm" ;
          exit (1) ;    
        }      
   }   
  edm::LogInfo ("IML") <<" [InvRingCalib][ctor] end of creator";
}

InvRingCalib::~InvRingCalib ( )

dtor

destructor

Definition at line 119 of file InvRingCalib.cc.

{
}

Member Function Documentation

void InvRingCalib::beginOfJob ( ) [virtual]

BeginOfJob.

Reimplemented from edm::EDLooperBase.

Definition at line 130 of file InvRingCalib.cc.

References isfirstcall_.

{
  isfirstcall_=true;

 
}

edm::EDLooper::Status InvRingCalib::duringLoop	(	const edm::Event &	iEvent,
		const edm::EventSetup &	iSetup
	)		`[virtual]`

duringLoop

Implements edm::EDLooper.

Definition at line 165 of file InvRingCalib.cc.

References DetId::det(), EBRegionDef(), DetId::Ecal, EcalBarrel, EcalEndcap, EERegionDef(), EERingDef(), HcalObjRepresent::Fill(), VFillMap::fillMap(), first, geometry, edm::EventSetup::get(), CaloMiscalibMapEcal::get(), edm::Event::getByLabel(), EcalCondObjectContainer< T >::getMap(), EcalClusterTools::getMaximum(), CaloGeometry::getValidDetIds(), EBDetId::ieta(), EBDetId::iphi(), isfirstcall_, edm::HandleBase::isValid(), EEDetId::ix(), EEDetId::iy(), edm::EDLooperBase::kContinue, m_barrelAlCa, m_barrelCells, m_barrelMap, m_EBcoeffFile, m_EEcoeffFile, m_ElectronLabel, m_endcapAlCa, m_endcapCells, m_endcapMap, m_IMACalibBlocks, m_MapFiller, m_mapFillerType, m_maxEnergyPerCrystal, m_minEnergyPerCrystal, m_recoWindowSidex, m_recoWindowSidey, m_RinginRegion, m_RingNumOfHits, m_xtalRegionId, m_xtalRing, siStripFEDMonitor_P5_cff::Max, dbtoconf::out, MiscalibReaderFromXML::parseXMLMiscalibFile(), CaloMiscalibMapEcal::prefillMap(), edm::Handle< T >::product(), DetId::rawId(), and EEDetId::zside().

{


   if (isfirstcall_){
    edm::LogInfo ("IML") << "[InvRingCalib][beginOfJob]" ;
    //gets the geometry from the event setup
    edm::ESHandle<CaloGeometry> geoHandle;
    iSetup.get<CaloGeometryRecord>().get(geoHandle);
    const CaloGeometry& geometry = *geoHandle;
    edm::LogInfo ("IML") <<"[InvRingCalib] Event Setup read";
    //fills a vector with all the cells
    m_barrelCells = geometry.getValidDetIds(DetId::Ecal, EcalBarrel);
    m_endcapCells = geometry.getValidDetIds(DetId::Ecal, EcalEndcap);
    //Defines the EB regions
    edm::LogInfo ("IML") <<"[InvRingCalib] Defining Barrel Regions";
    EBRegionDef();
    //Defines what is a ring in the EE
    edm::LogInfo ("IML") <<"[InvRingCalib] Defining endcap Rings";
    EERingDef(iSetup);
    //Defines the regions in the EE
    edm::LogInfo ("IML") <<"[InvRingCalib] Defining endcap Regions";
    EERegionDef();
    if (m_mapFillerType == "Cluster") m_MapFiller= new ClusterFillMap (
                                                                       m_recoWindowSidex ,m_recoWindowSidey ,
                                                                       m_xtalRegionId ,m_minEnergyPerCrystal ,
                                                                       m_maxEnergyPerCrystal , m_RinginRegion ,
                                                                       & m_barrelMap ,
                                                                       & m_endcapMap ); 
    if (m_mapFillerType == "Matrix") m_MapFiller = new MatrixFillMap (
                                                                      m_recoWindowSidex ,m_recoWindowSidey ,
                                                                      m_xtalRegionId , m_minEnergyPerCrystal ,
                                                                      m_maxEnergyPerCrystal , m_RinginRegion ,
                                                                      & m_barrelMap ,
                                                                      & m_endcapMap); 
    edm::LogInfo ("IML") <<"[InvRingCalib] Initializing the coeffs";
    //Sets the initial coefficients to 1.
    //Graphs to check ring, regions and so on, not needed in the final version
    TH2F EBRegion ("EBRegion","EBRegion",171,-85,86,360,1,361);
    TH2F EBRing ("EBRing","EBRing",171,-85,86,360,1,361);
    for (std::vector<DetId>::const_iterator it= m_barrelCells.begin();
         it!= m_barrelCells.end(); 
         ++it )
      {
        EBDetId eb (*it);
        EBRing.Fill(eb.ieta(),eb.iphi(),m_RinginRegion[it->rawId()]);
        EBRegion.Fill(eb.ieta(),eb.iphi(),m_xtalRegionId[it->rawId()]);
      }

    TH2F EEPRegion ("EEPRegion", "EEPRegion",100,1,101,100,1,101);
    TH2F EEPRing ("EEPRing", "EEPRing",100,1,101,100,1,101);
    TH2F EEPRingReg ("EEPRingReg", "EEPRingReg",100,1,101,100,1,101);
    TH2F EEMRegion ("EEMRegion", "EEMRegion",100,1,101,100,1,101);
    TH2F EEMRing ("EEMRing", "EEMRing",100,1,101,100,1,101);
    TH2F EEMRingReg ("EEMRingReg", "EEMRingReg",100,1,101,100,1,101);
    //  TH1F eta ("eta","eta",250,-85,165);
    for (std::vector<DetId>::const_iterator it = m_endcapCells.begin();
         it!= m_endcapCells.end();
         ++it)
      {
        EEDetId ee (*it);
        if (ee.zside()>0)
          {
            EEPRegion.Fill(ee.ix(),ee.iy(),m_xtalRegionId[ee.rawId()]);
            EEPRing.Fill(ee.ix(),ee.iy(),m_xtalRing[ee.rawId()]);
            EEPRingReg.Fill(ee.ix(),ee.iy(),m_RinginRegion[ee.rawId()]);
          }
        if (ee.zside()<0)
          {
            EEMRegion.Fill(ee.ix(),ee.iy(),m_xtalRegionId[ee.rawId()]);
            EEMRing.Fill(ee.ix(),ee.iy(),m_xtalRing[ee.rawId()]);
            EEMRingReg.Fill(ee.ix(),ee.iy(),m_RinginRegion[ee.rawId()]);
          }    
      } 

    //  for (std::map<int,float>::iterator it=m_eta.begin();
    //        it!=m_eta.end();++it)
    //             eta.Fill(it->first,it->second);
    TFile out ("EBZone.root", "recreate");
    EBRegion.Write();
    EBRing.Write();
    EEPRegion.Write();
    EEPRing.Write();
    EEPRingReg.Write();
    EEMRegion.Write();
    EEMRing.Write();
    //  eta.Write();
    EEMRingReg.Write();
    out.Close();
    edm::LogInfo ("IML") <<"[InvRingCalib] Start to acquire the coeffs";
    CaloMiscalibMapEcal EBmap;
    EBmap.prefillMap ();
    MiscalibReaderFromXMLEcalBarrel barrelreader (EBmap);
    if (!m_EBcoeffFile.empty()) barrelreader.parseXMLMiscalibFile (m_EBcoeffFile);
    EcalIntercalibConstants costants (EBmap.get());
    m_barrelMap = costants.getMap();
    CaloMiscalibMapEcal EEmap ;   
    EEmap.prefillMap ();
    MiscalibReaderFromXMLEcalEndcap endcapreader (EEmap);
    if (!m_EEcoeffFile.empty()) endcapreader.parseXMLMiscalibFile (m_EEcoeffFile) ;
    EcalIntercalibConstants EEcostants (EEmap.get());
    m_endcapMap = EEcostants.getMap();

    isfirstcall_=false;
  } // if isfirstcall






  //gets the barrel recHits
  double pSubtract = 0.;
  double pTk = 0.;
  const EcalRecHitCollection* barrelHitsCollection = 0;
  edm::Handle<EBRecHitCollection> barrelRecHitsHandle ;
  iEvent.getByLabel (m_barrelAlCa, barrelRecHitsHandle) ;
  barrelHitsCollection = barrelRecHitsHandle.product () ;

 if (!barrelRecHitsHandle.isValid ()) {
     edm::LogError ("IML") << "[EcalEleCalibLooper] barrel rec hits not found" ;
     return  kContinue ;
    }
  //gets the endcap recHits
  const EcalRecHitCollection* endcapHitsCollection = 0;
  edm::Handle<EERecHitCollection> endcapRecHitsHandle ;
  iEvent.getByLabel (m_endcapAlCa, endcapRecHitsHandle) ;
  endcapHitsCollection = endcapRecHitsHandle.product () ;

 if (!endcapRecHitsHandle.isValid ()) {  
     edm::LogError ("IML") << "[EcalEleCalibLooper] endcap rec hits not found" ; 
     return kContinue;
   }

  //gets the electrons
  edm::Handle<reco::GsfElectronCollection> pElectrons;
  iEvent.getByLabel(m_ElectronLabel,pElectrons);

 if (!pElectrons.isValid ()) {
     edm::LogError ("IML")<< "[EcalEleCalibLooper] electrons not found" ;
     return kContinue;
   }

  //loops over the electrons in the event
  for (reco::GsfElectronCollection::const_iterator eleIt = pElectrons->begin();
       eleIt != pElectrons->end();
       ++eleIt )
    {
      pSubtract =0;
      pTk=eleIt->trackMomentumAtVtx().R();
      std::map<int , double> xtlMap;
      DetId Max=0; 
      if (fabs(eleIt->eta()<1.49))
             Max = EcalClusterTools::getMaximum(eleIt->superCluster()->hitsAndFractions(),barrelHitsCollection).first;
      else 
             Max = EcalClusterTools::getMaximum(eleIt->superCluster()->hitsAndFractions(),endcapHitsCollection).first;
      if (Max.det()==0) continue;
       m_MapFiller->fillMap(eleIt->superCluster ()->hitsAndFractions (),Max, 
                           barrelHitsCollection,endcapHitsCollection, xtlMap,pSubtract);
      if (m_xtalRegionId[Max.rawId()]==-1) continue;
      pSubtract += eleIt->superCluster()->preshowerEnergy() ;
      ++m_RingNumOfHits[m_xtalRing[Max.rawId()]];
      //fills the calibBlock 
      m_IMACalibBlocks.at(m_xtalRegionId[Max.rawId()])->Fill (
          xtlMap.begin(), xtlMap.end(),pTk,pSubtract
        ) ;
    }
  return  kContinue;
} //end of duringLoop

int InvRingCalib::EBRegId ( const int ieta ) [private]

Gives back in which region you are:

gives the region Id given ieta

Definition at line 609 of file InvRingCalib.cc.

References m_etaEnd, and m_Reg.

Referenced by EBRegionDef().

{
 if (ieta<m_etaStart || ieta>=m_etaEnd) return -1;
 else return (m_Reg[ieta]);
}

void InvRingCalib::EBRegionDef ( ) [private]

Defines the regions in the barrel.

Definition at line 620 of file InvRingCalib.cc.

References EBRegId(), m_barrelCells, m_etaStart, m_etaWidth, m_RinginRegion, m_xtalRegionId, m_xtalRing, and RegPrepare().

Referenced by duringLoop().

{
  RegPrepare();
  for (std::vector<DetId>::const_iterator it=m_barrelCells.begin();
        it!=m_barrelCells.end();++it)
  {
    EBDetId eb (it->rawId());
    m_xtalRing[eb.rawId()] = eb.ieta() ;
    m_xtalRegionId[eb.rawId()] = EBRegId (eb.ieta()); 
    if (m_xtalRegionId[eb.rawId()]==-1) continue;
    m_RinginRegion[eb.rawId()] = (eb.ieta() - m_etaStart)% m_etaWidth; 
  }
}

int InvRingCalib::EBRegionNum ( ) const [private]

Number of regions in EB.

number of Ring in EB

Definition at line 586 of file InvRingCalib.cc.

References m_etaEnd, m_etaStart, and m_etaWidth.

Referenced by EERegionDef(), and InvRingCalib().

{
  if ((m_etaEnd*m_etaStart)>0)
   return ((m_etaEnd - m_etaStart )/m_etaWidth); 
  
  if ((m_etaEnd*m_etaStart)<0)
   return ((m_etaEnd - m_etaStart-1 )/m_etaWidth); 
  
  return 0;
}

int InvRingCalib::EERegId ( int id ) [private]

gives back in which region of the endcap you are.

Gives the Id of the region given the id of the xtal.

Definition at line 518 of file InvRingCalib.cc.

References m_EEZone, m_endRing, m_etaWidth, m_RinginRegion, m_startRing, m_xtalRing, relativeConstraints::ring, and EEDetId::zside().

Referenced by EERegionDef().

{
   int reg;
   int ring;
   EEDetId ee (id);
  //sets the reg to -1 if the ring doesn't exist or is outside the region of interest 
   if (m_xtalRing[id] == -1) return -1;
  //avoid the calibration in the wrong zside
   if (m_EEZone == 1 ){
   if (ee.zside()<0) return -1;
   ring = m_xtalRing[id]-86;
   if(ring >=m_endRing) return -1;
   if (ring<m_startRing) return -1;
   reg = (ring -m_startRing) / m_etaWidth;
   m_RinginRegion[id]=(ring -m_startRing) % m_etaWidth;
   return reg;
   }
   if (m_EEZone == -1){
   if (ee.zside()>0) return -1;
   ring = m_xtalRing[id] -125;
   if(ring >=m_endRing) return -1;
   if (ring<m_startRing) return -1;
   reg = (ring -m_startRing) / m_etaWidth;
   m_RinginRegion[id]=(ring -m_startRing) % m_etaWidth;
   return reg;
   }
   if (ee.zside()>0) ring=m_xtalRing[id]-86;
     else ring = m_xtalRing[id]-125;
   if(ring >=m_endRing) return -1;
   if (ring<m_startRing) return -1;
   reg = (ring -m_startRing) / m_etaWidth;
   m_RinginRegion[id]=(ring -m_startRing) % m_etaWidth;
   return reg;
}

void InvRingCalib::EERegionDef ( ) [private]

Defines the regions in the endcap.

Loops over all the endcap xtals and sets for each xtal the value of the region the xtal is in, and the ringNumber inside the region

Definition at line 555 of file InvRingCalib.cc.

References EBRegionNum(), EERegId(), EERegionNum(), m_endcapCells, m_xtalRegionId, and EEDetId::zside().

Referenced by duringLoop().

{
int reg;
for (std::vector<DetId>::const_iterator endcapIt=m_endcapCells.begin();
     endcapIt!=m_endcapCells.end();++endcapIt){
      EEDetId ee(*endcapIt);
      reg = EERegId(endcapIt->rawId());
      //If the ring is not of interest saves only the region Id(-1)
      if(reg==-1) 
         m_xtalRegionId[endcapIt->rawId()]=reg;
      //sums the number of region in EB or EB+EE to have different regionsId in different regions 
      else {
      if (ee.zside()>0)reg += EBRegionNum();
      else reg += EBRegionNum()+EERegionNum();
      m_xtalRegionId[endcapIt->rawId()]=reg;
   }
  }
}

int InvRingCalib::EERegionNum ( ) const [inline, private]

The number of regions in EE.

Number of Regions in EE.

Definition at line 579 of file InvRingCalib.cc.

References m_endRing, m_etaWidth, and m_startRing.

Referenced by EERegionDef(), and InvRingCalib().

{
  return ((m_endRing - m_startRing)/m_etaWidth);
}

void InvRingCalib::EERingDef ( const edm::EventSetup & iSetup ) [private]

Defines the rins in the endcap.

EE ring definition.

Definition at line 455 of file InvRingCalib.cc.

References DetId::Ecal, EcalEndcap, geometry, edm::EventSetup::get(), CaloSubdetectorGeometry::getGeometry(), CaloCellGeometry::getPosition(), CaloGeometry::getSubdetectorGeometry(), m_cellPhi, m_cellPos, m_endcapCells, m_xtalRing, PV3DBase< T, PVType, FrameType >::phi(), DetId::rawId(), relativeConstraints::ring, EEDetId::validDetId(), EEDetId::XYMODE, and EEDetId::zside().

Referenced by duringLoop().

{
 //Gets the Handle for the geometry from the eventSetup
 edm::ESHandle<CaloGeometry> geoHandle;
 iSetup.get<CaloGeometryRecord>().get(geoHandle);
 //Gets the geometry of the endcap
 const CaloGeometry& geometry = *geoHandle;
 const CaloSubdetectorGeometry *endcapGeometry = geometry.getSubdetectorGeometry(DetId::Ecal, EcalEndcap);
// const CaloSubdetectorGeometry *barrelGeometry = geometry.getSubdetectorGeometry(DetId::Ecal, EcalBarrel);
 //for every xtal gets the position Vector and the phi position
 
// for (std::vector<DetId>::const_iterator barrelIt = m_barrelCells.begin();
//    barrelIt!=m_barrelCells.end();
//    ++barrelIt) {
//     const CaloCellGeometry *cellGeometry = barrelGeometry->getGeometry(*barrelIt);
//     GlobalPoint point;
//     EBDetId eb (*barrelIt);
//     point=cellGeometry->getPosition();
//     m_eta[eb.ieta()]=point.eta() ;    //cellGeometry->getPosition().eta();
//    }
 for (std::vector<DetId>::const_iterator endcapIt = m_endcapCells.begin();
    endcapIt!=m_endcapCells.end();
    ++endcapIt) {
     const CaloCellGeometry *cellGeometry = endcapGeometry->getGeometry(*endcapIt);
     m_cellPos[endcapIt->rawId()] = cellGeometry->getPosition();
     m_cellPhi[endcapIt->rawId()] = cellGeometry->getPosition().phi();
  }
 //takes the first 39 xtals at a fixed y varying the x coordinate and saves their eta coordinate 
 float eta_ring[39];
 for (int ring=0; ring<39; ring++) 
        if (EEDetId::validDetId(ring,50,1)){  
          EEDetId det = EEDetId (ring,50,1,EEDetId::XYMODE);
          eta_ring[ring]=m_cellPos[det.rawId()].eta();
          }
 //defines the bonduary of the rings as the average eta of a xtal
 double etaBonduary[40];
 etaBonduary[0]=1.49;
 etaBonduary[39]=4.0;
 for (int ring=1; ring<39; ++ring)
       etaBonduary[ring]=(eta_ring[ring]+eta_ring[ring-1])/2.;
 //assign to each xtal a ring number
 int CRing;
 for (int ring=0; ring<39; ring++) 
   for (std::vector<DetId>::const_iterator endcapIt=m_endcapCells.begin();
        endcapIt!=m_endcapCells.end();++endcapIt){
     if (fabs(m_cellPos[endcapIt->rawId()].eta())>etaBonduary[ring] &&
         fabs(m_cellPos[endcapIt->rawId()].eta())<etaBonduary[ring+1])
          {
              EEDetId ee(*endcapIt);
              if (ee.zside()>0) CRing=ring + 86; 
              else CRing = ring + 125;
              m_xtalRing[endcapIt->rawId()]=CRing;
//              m_eta[CRing]=m_cellPos[endcapIt->rawId()].eta();
          }    
      }
 return;
}

void InvRingCalib::endOfJob ( ) [virtual]

Reimplemented from edm::EDLooperBase.

Definition at line 426 of file InvRingCalib.cc.

References EcalBarrel, EcalEndcap, m_barrelCells, m_barrelMap, m_endcapCells, m_endcapMap, and calibXMLwriter::writeLine().

{

 edm::LogInfo ("IML") << "[InvMatrixCalibLooper][endOfJob] saving calib coeffs" ;
 calibXMLwriter barrelWriter(EcalBarrel);
 calibXMLwriter endcapWriter(EcalEndcap);
 for (std::vector<DetId>::const_iterator barrelIt =m_barrelCells.begin(); 
       barrelIt!=m_barrelCells.end(); 
       ++barrelIt) {
            EBDetId eb (*barrelIt);
            barrelWriter.writeLine(eb,m_barrelMap[eb]);
          }
 for (std::vector<DetId>::const_iterator endcapIt = m_endcapCells.begin();
     endcapIt!=m_endcapCells.end();
     ++endcapIt) {
          EEDetId ee (*endcapIt);
          endcapWriter.writeLine(ee,m_endcapMap[ee]);
        }
}

edm::EDLooper::Status InvRingCalib::endOfLoop	(	const edm::EventSetup &	,
		unsigned int	iCounter
	)		`[virtual]`

Called after the system has finished one loop over the events. Thar argument is a count of how many loops have been processed before this loo. For the first time through the events the argument will be 0.

Implements edm::EDLooperBase.

Definition at line 342 of file InvRingCalib.cc.

References lut2db_cfg::filename, edm::EDLooperBase::kContinue, edm::EDLooperBase::kStop, m_barrelCells, m_barrelMap, m_endcapCells, m_endcapMap, m_IMACalibBlocks, m_loops, m_maxCoeff, m_minCoeff, m_RinginRegion, m_RingNumOfHits, m_usingBlockSolver, m_xtalRegionId, m_xtalRing, and dbtoconf::out.

{
   std::map<int,double> InterRings;
  edm::LogInfo ("IML") << "[InvMatrixCalibLooper][endOfLoop] Start to invert the matrixes" ;
  //call the autoexplaining "solve" method for every calibBlock
  for (std::vector<VEcalCalibBlock *>::iterator calibBlock=m_IMACalibBlocks.begin();
       calibBlock!=m_IMACalibBlocks.end();
       ++calibBlock)
    (*calibBlock)->solve(m_usingBlockSolver,m_minCoeff,m_maxCoeff);

  edm::LogInfo("IML") << "[InvRingLooper][endOfLoop] Starting to write the coeffs";
  TH1F *coeffDistr = new TH1F("coeffdistr","coeffdistr",100 ,0.7,1.4);
  TH1F *coeffMap = new TH1F("coeffRingMap","coeffRingMap",250,-85,165);
  TH1F *ringDistr = new TH1F("ringDistr","ringDistr",250,-85,165);
  TH1F *RingFill = new TH1F("RingFill","RingFill",250,-85,165);
  for(std::map<int,int>::const_iterator it=m_xtalRing.begin();
      it!=m_xtalRing.end();
      ++it)
    ringDistr->Fill(it->second+0.1);

  int ID;
  std::map<int,int> flag;
  for(std::map<int,int>::const_iterator it=m_xtalRing.begin();
      it!=m_xtalRing.end();
      ++it)
    flag[it->second]=0;
  
  for (std::vector<DetId>::const_iterator it=m_barrelCells.begin();
       it!=m_barrelCells.end();
       ++it)
    { 
      ID= it->rawId();
      if (m_xtalRegionId[ID]==-1) continue;
      if (flag[m_xtalRing[ID]]) continue;
      flag[m_xtalRing[ID]] =1;
      RingFill->Fill(m_xtalRing[ID],m_RingNumOfHits[m_xtalRing[ID]]);
      InterRings[m_xtalRing[ID]] = m_IMACalibBlocks.at(m_xtalRegionId[ID])->at(m_RinginRegion[ID]);
      coeffMap->Fill (m_xtalRing[ID]+0.1,InterRings[m_xtalRing[ID]]);
      coeffDistr->Fill(InterRings[m_xtalRing[ID]]);
    }

  for (std::vector<DetId>::const_iterator it=m_endcapCells.begin();
       it!=m_endcapCells.end();
       ++it)
    { 
      ID= it->rawId();
      if (m_xtalRegionId[ID]==-1) continue;
      if (flag[m_xtalRing[ID]]) continue;
      flag[m_xtalRing[ID]]= 1;
      InterRings[m_xtalRing[ID]] = m_IMACalibBlocks.at(m_xtalRegionId[ID])->at(m_RinginRegion[ID]);
      RingFill->Fill(m_xtalRing[ID],m_RingNumOfHits[m_xtalRing[ID]]);
      coeffMap->Fill (m_xtalRing[ID],InterRings[m_xtalRing[ID]]);
      coeffDistr->Fill(InterRings[m_xtalRing[ID]]);
                
    } 
   
  char filename[80];
  sprintf(filename,"coeff%d.root",iCounter);
  TFile out(filename,"recreate");    
  coeffDistr->Write();
  coeffMap->Write();
  ringDistr->Write();
  RingFill->Write();
  out.Close();
  for (std::vector<DetId>::const_iterator it=m_barrelCells.begin();
       it!=m_barrelCells.end();
       ++it){
         m_barrelMap[*it]*=InterRings[m_xtalRing[it->rawId()]];
  }
  for (std::vector<DetId>::const_iterator it=m_endcapCells.begin();
       it!=m_endcapCells.end();
       ++it)
          m_endcapMap[*it]*=InterRings[m_xtalRing[it->rawId()]];
  if (iCounter < m_loops-1 ) return kContinue ;
  else return kStop; 
}

void InvRingCalib::RegPrepare ( ) [private]

Prepares the EB regions;.

Divides the barrel in region, necessary to take into account the missing 0 xtal

Definition at line 598 of file InvRingCalib.cc.

References i, gen::k, m_etaEnd, m_etaStart, m_etaWidth, and m_Reg.

Referenced by EBRegionDef().

{
 int k=0;
 for (int i = m_etaStart;i<m_etaEnd;++i)
 {
  if (i==0) continue;
  m_Reg[i]=k/m_etaWidth;
  ++k;
 }
}

void InvRingCalib::startingNewLoop ( unsigned int ciclo ) [virtual]

startingNewLoop

Implements edm::EDLooperBase.

Definition at line 142 of file InvRingCalib.cc.

References m_IMACalibBlocks, m_RingNumOfHits, m_xtalRing, and relativeConstraints::ring.

{
    edm::LogInfo ("IML") << "[InvMatrixCalibLooper][Start] entering loop " << ciclo;
    for (std::vector<VEcalCalibBlock *>::iterator calibBlock = m_IMACalibBlocks.begin () ;
         calibBlock != m_IMACalibBlocks.end () ;
         ++calibBlock)
      {
        //LP empties the energies vector, to fill DuringLoop.
        (*calibBlock)->reset () ;
      }
   for (std::map<int,int>::const_iterator ring=m_xtalRing.begin();
        ring!=m_xtalRing.end();
        ++ring)
           m_RingNumOfHits[ring->second]=0;
   return ;
}