composite_calorimeter/src/CCaloOrganization.cc

0001 //
0002 // ********************************************************************
0003 // * License and Disclaimer                                           *
0004 // *                                                                  *
0005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
0006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
0007 // * conditions of the Geant4 Software License,  included in the file *
0008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
0009 // * include a list of copyright holders.                             *
0010 // *                                                                  *
0011 // * Neither the authors of this software system, nor their employing *
0012 // * institutes,nor the agencies providing financial support for this *
0013 // * work  make  any representation or  warranty, express or implied, *
0014 // * regarding  this  software system or assume any liability for its *
0015 // * use.  Please see the license in the file  LICENSE  and URL above *
0016 // * for the full disclaimer and the limitation of liability.         *
0017 // *                                                                  *
0018 // * This  code  implementation is the result of  the  scientific and *
0019 // * technical work of the GEANT4 collaboration.                      *
0020 // * By using,  copying,  modifying or  distributing the software (or *
0021 // * any work based  on the software)  you  agree  to acknowledge its *
0022 // * use  in  resulting  scientific  publications,  and indicate your *
0023 // * acceptance of all terms of the Geant4 Software license.          *
0024 // ********************************************************************
0025 //
0026 ///////////////////////////////////////////////////////////////////////////////
0027 // File: CCaloOrganization.cc
0028 // Description: Packing, unpacking and other related utilities for
0029 //              calorimetric numbering schema
0030 ///////////////////////////////////////////////////////////////////////////////
0031 #include "CCaloOrganization.hh"
0032
0033 //#define debug
0034
0035 unsigned int CCaloOrganization::packindex(G4int det, G4int z, G4int eta,
0036                                           G4int phi) const {
0037   //So this is the actual encoding of the index:
0038   //top 4 bits encode Detector type
0039   //
0040   // Should work for all calorimeter with no depth information
0041
0042   unsigned int idx=(det&15)<<28;   //bits 28-31   (21-23 are free for now)
0043   idx+=(((z+1)/2)&1)<<20;          //bits 20
0044   idx+=(eta&1023)<<10;             //bits 10-19
0045   idx+=(phi&1023);                 //bits  0-9
0046 #ifdef debug
0047   G4cout << " ECAL packing " << det << " " << z << " " << eta << " " << phi
0048          << "  into " << idx << G4endl;
0049 #endif
0050   return idx;
0051 }
0052
0053 unsigned int CCaloOrganization::packindex(G4int det, G4int depth, G4int z, G4int eta,
0054                                           G4int phi) const {
0055   //So this is the actual encoding of the index:
0056   //top 4 bits encode Detector type
0057   //next 4 bits encode depth information
0058   // Should work for all calorimeter with no depth information
0059
0060   unsigned int idx=(det&15)<<28;  //bits 28-31   (21-23 are free for now)
0061   idx+=(depth&15)<<24;            //bits 24-27
0062   idx+=(z&1)<<20;                 //bits 20
0063   idx+=(eta&1023)<<10;            //bits 10-19
0064   idx+=(phi&1023);                //bits  0-9
0065 #ifdef debug
0066   G4cout << " HCAL packing " << det << " " << depth << " " << z << " " << eta
0067          << " " << phi  << "  into " << idx << G4endl;
0068 #endif
0069   return idx;
0070 }
0071
0072
0073 void CCaloOrganization::unpackindex(const unsigned int& idx, G4int& det, G4int& z,
0074                                     G4int& eta, G4int& phi) const {
0075   det = (idx>>28)&15;
0076   z   = (idx>>20)&1;
0077   z   = 2*z-1;
0078   eta = (idx>>10)&1023;
0079   phi = (idx&1023);
0080 }
0081
0082
0083 void CCaloOrganization::unpackindex(const unsigned int& idx, G4int& det,
0084                                     G4int& depth, G4int& z, G4int& eta,
0085                                     G4int& phi) const {
0086   det = (idx>>28)&15;
0087   depth=(idx>>24)&15;
0088   z   = (idx>>20)&1;
0089   eta = (idx>>10)&1023;
0090   phi = (idx&1023);
0091 }
0092
0093
0094 G4int CCaloOrganization::getUnitWithMaxEnergy(std::map<G4int,G4float,std::less<G4int> >& themap){
0095
0096   //look for max
0097   G4int UnitWithMaxEnergy = 0;
0098   G4float maxEnergy = 0.;
0099
0100   for(std::map<G4int,G4float,std::less<G4int> >::iterator iter = themap.begin();
0101       iter != themap.end(); iter++){
0102
0103     if(        maxEnergy < (*iter).second) {
0104       maxEnergy = (*iter).second;
0105       UnitWithMaxEnergy = (*iter).first;
0106     }
0107   }
0108 #ifdef debug
0109   G4cout << " *** max energy of " << maxEnergy << " MeV was found in Unit id "
0110          << UnitWithMaxEnergy;
0111   G4int det,z,eta,phi;
0112   unpackindex(UnitWithMaxEnergy, det, z, eta, phi);
0113   G4cout << " corresponding to z= " << z << " eta= " << eta << " phi = " << phi
0114          << G4endl;
0115 #endif
0116   return UnitWithMaxEnergy;
0117
0118 }
0119
0120
0121 G4float CCaloOrganization::energyInMatrix(G4int nCellInEta, G4int nCellInPhi,
0122                                         G4int crystalWithMaxEnergy,
0123                                         std::map<G4int,G4float,std::less<G4int> >& themap){
0124
0125   G4int det,z,eta,phi;
0126   this->unpackindex(crystalWithMaxEnergy, det, z, eta, phi);
0127 #ifdef debug
0128   G4int ncristals=0;
0129 #endif
0130   G4int goBackInEta = nCellInEta/2;
0131   G4int goBackInPhi = nCellInPhi/2;
0132   G4int startEta = eta-goBackInEta;
0133   G4int startPhi = phi-goBackInPhi;
0134
0135   G4float totalEnergy = 0.;
0136
0137   for(G4int ieta=startEta; ieta<startEta+nCellInEta; ieta++){
0138     for(G4int iphi=startPhi; iphi<startPhi+nCellInPhi; iphi++){
0139
0140       G4int index = this->packindex(det,z,ieta,iphi);
0141       totalEnergy += themap[index];
0142 #ifdef debug
0143       ncristals+=1;
0144       G4cout << "ieta - iphi - E = " << ieta << "  " << iphi << " "
0145              << themap[index] << G4endl;
0146 #endif
0147     }
0148   }
0149
0150 #ifdef debug
0151   G4cout << "Energy in " << nCellInEta << " cells in eta times "
0152          << nCellInPhi << " cells in phi matrix = " << totalEnergy
0153          << " for " << ncristals << " crystals" << G4endl;
0154 #endif
0155   return totalEnergy;
0156
0157 }