| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer << 3 // * DISCLAIMER *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th << 5 // * The following disclaimer summarizes all the specific disclaimers *
6 // * the Geant4 Collaboration. It is provided << 6 // * of contributors to this software. The specific disclaimers,which *
7 // * conditions of the Geant4 Software License << 7 // * govern, are listed with their locations in: *
8 // * LICENSE and available at http://cern.ch/ << 8 // * http://cern.ch/geant4/license *
9 // * include a list of copyright holders. <<
10 // * 9 // * *
11 // * Neither the authors of this software syst 10 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 11 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 12 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 13 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file << 14 // * use. *
16 // * for the full disclaimer and the limitatio <<
17 // * 15 // * *
18 // * This code implementation is the result << 16 // * This code implementation is the intellectual property of the *
19 // * technical work of the GEANT4 collaboratio << 17 // * GEANT4 collaboration. *
20 // * By using, copying, modifying or distri << 18 // * By copying, distributing or modifying the Program (or any work *
21 // * any work based on the software) you ag << 19 // * based on the Program) you indicate your acceptance of this *
22 // * use in resulting scientific publicati << 20 // * statement, and all its terms. *
23 // * acceptance of all terms of the Geant4 Sof <<
24 // ******************************************* 21 // ********************************************************************
25 // 22 //
26 ////////////////////////////////////////////// 23 ///////////////////////////////////////////////////////////////////////////////
27 // File: CCaloOrganization.cc 24 // File: CCaloOrganization.cc
28 // Description: Packing, unpacking and other r 25 // Description: Packing, unpacking and other related utilities for
29 // calorimetric numbering schema 26 // calorimetric numbering schema
30 ////////////////////////////////////////////// 27 ///////////////////////////////////////////////////////////////////////////////
31 #include "CCaloOrganization.hh" 28 #include "CCaloOrganization.hh"
32 29
33 //#define debug 30 //#define debug
34 31
35 unsigned int CCaloOrganization::packindex(G4in << 32 unsigned int CCaloOrganization::packindex(int det, int z, int eta,
36 G4in << 33 int phi) const {
37 //So this is the actual encoding of the inde 34 //So this is the actual encoding of the index:
38 //top 4 bits encode Detector type 35 //top 4 bits encode Detector type
39 // 36 //
40 // Should work for all calorimeter with no d 37 // Should work for all calorimeter with no depth information
41 38
42 unsigned int idx=(det&15)<<28; //bits 28-3 39 unsigned int idx=(det&15)<<28; //bits 28-31 (21-23 are free for now)
43 idx+=(((z+1)/2)&1)<<20; //bits 20 40 idx+=(((z+1)/2)&1)<<20; //bits 20
44 idx+=(eta&1023)<<10; //bits 10-1 41 idx+=(eta&1023)<<10; //bits 10-19
45 idx+=(phi&1023); //bits 0-9 42 idx+=(phi&1023); //bits 0-9
46 #ifdef debug 43 #ifdef debug
47 G4cout << " ECAL packing " << det << " " << 44 G4cout << " ECAL packing " << det << " " << z << " " << eta << " " << phi
48 << " into " << idx << G4endl; << 45 << " into " << idx << G4endl;
49 #endif 46 #endif
50 return idx; 47 return idx;
51 } 48 }
52 49
53 unsigned int CCaloOrganization::packindex(G4in << 50 unsigned int CCaloOrganization::packindex(int det, int depth, int z, int eta,
54 G4in << 51 int phi) const {
55 //So this is the actual encoding of the inde 52 //So this is the actual encoding of the index:
56 //top 4 bits encode Detector type 53 //top 4 bits encode Detector type
57 //next 4 bits encode depth information 54 //next 4 bits encode depth information
58 // Should work for all calorimeter with no d 55 // Should work for all calorimeter with no depth information
59 56
60 unsigned int idx=(det&15)<<28; //bits 28-31 57 unsigned int idx=(det&15)<<28; //bits 28-31 (21-23 are free for now)
61 idx+=(depth&15)<<24; //bits 24-27 58 idx+=(depth&15)<<24; //bits 24-27
62 idx+=(z&1)<<20; //bits 20 59 idx+=(z&1)<<20; //bits 20
63 idx+=(eta&1023)<<10; //bits 10-19 60 idx+=(eta&1023)<<10; //bits 10-19
64 idx+=(phi&1023); //bits 0-9 61 idx+=(phi&1023); //bits 0-9
65 #ifdef debug 62 #ifdef debug
66 G4cout << " HCAL packing " << det << " " << 63 G4cout << " HCAL packing " << det << " " << depth << " " << z << " " << eta
67 << " " << phi << " into " << idx << << 64 << " " << phi << " into " << idx << G4endl;
68 #endif 65 #endif
69 return idx; 66 return idx;
70 } 67 }
71 68
72 69
73 void CCaloOrganization::unpackindex(const unsi << 70 void CCaloOrganization::unpackindex(const unsigned int& idx, int& det, int& z,
74 G4int& eta << 71 int& eta, int& phi) const {
75 det = (idx>>28)&15; 72 det = (idx>>28)&15;
76 z = (idx>>20)&1; 73 z = (idx>>20)&1;
77 z = 2*z-1; 74 z = 2*z-1;
78 eta = (idx>>10)&1023; 75 eta = (idx>>10)&1023;
79 phi = (idx&1023); 76 phi = (idx&1023);
80 } 77 }
81 78
82 79
83 void CCaloOrganization::unpackindex(const unsi << 80 void CCaloOrganization::unpackindex(const unsigned int& idx, int& det,
84 G4int& dep << 81 int& depth, int& z, int& eta,
85 G4int& phi << 82 int& phi) const {
86 det = (idx>>28)&15; 83 det = (idx>>28)&15;
87 depth=(idx>>24)&15; 84 depth=(idx>>24)&15;
88 z = (idx>>20)&1; 85 z = (idx>>20)&1;
89 eta = (idx>>10)&1023; 86 eta = (idx>>10)&1023;
90 phi = (idx&1023); 87 phi = (idx&1023);
91 } 88 }
92 89
93 90
94 G4int CCaloOrganization::getUnitWithMaxEnergy( << 91 int CCaloOrganization::getUnitWithMaxEnergy(std::map<int,float,std::less<int> >& themap){
95 92
96 //look for max 93 //look for max
97 G4int UnitWithMaxEnergy = 0; << 94 int UnitWithMaxEnergy = 0;
98 G4float maxEnergy = 0.; << 95 float maxEnergy = 0.;
99 << 96
100 for(std::map<G4int,G4float,std::less<G4int> << 97 for(std::map<int,float,std::less<int> >::iterator iter = themap.begin();
101 iter != themap.end(); iter++){ 98 iter != themap.end(); iter++){
102 << 99
103 if( maxEnergy < (*iter).second) { << 100 if( maxEnergy < (*iter).second) {
104 maxEnergy = (*iter).second; << 101 maxEnergy = (*iter).second;
105 UnitWithMaxEnergy = (*iter).first; 102 UnitWithMaxEnergy = (*iter).first;
106 } << 103 }
107 } << 104 }
108 #ifdef debug 105 #ifdef debug
109 G4cout << " *** max energy of " << maxEnergy 106 G4cout << " *** max energy of " << maxEnergy << " MeV was found in Unit id "
110 << UnitWithMaxEnergy; << 107 << UnitWithMaxEnergy;
111 G4int det,z,eta,phi; << 108 int det,z,eta,phi;
112 unpackindex(UnitWithMaxEnergy, det, z, eta, 109 unpackindex(UnitWithMaxEnergy, det, z, eta, phi);
113 G4cout << " corresponding to z= " << z << " 110 G4cout << " corresponding to z= " << z << " eta= " << eta << " phi = " << phi
114 << G4endl; << 111 << G4endl;
115 #endif 112 #endif
116 return UnitWithMaxEnergy; 113 return UnitWithMaxEnergy;
117 114
118 } 115 }
119 116
120 117
121 G4float CCaloOrganization::energyInMatrix(G4in << 118 float CCaloOrganization::energyInMatrix(int nCellInEta, int nCellInPhi,
122 G4int << 119 int crystalWithMaxEnergy,
123 std::m << 120 std::map<int,float,std::less<int> >& themap){
124 121
125 G4int det,z,eta,phi; << 122 int det,z,eta,phi;
126 this->unpackindex(crystalWithMaxEnergy, det, 123 this->unpackindex(crystalWithMaxEnergy, det, z, eta, phi);
127 #ifdef debug << 124 int ncristals=0;
128 G4int ncristals=0; << 125
129 #endif << 126 int goBackInEta = nCellInEta/2;
130 G4int goBackInEta = nCellInEta/2; << 127 int goBackInPhi = nCellInPhi/2;
131 G4int goBackInPhi = nCellInPhi/2; << 128 int startEta = eta-goBackInEta;
132 G4int startEta = eta-goBackInEta; << 129 int startPhi = phi-goBackInPhi;
133 G4int startPhi = phi-goBackInPhi; <<
134 130
135 G4float totalEnergy = 0.; << 131 float totalEnergy = 0.;
136 132
137 for(G4int ieta=startEta; ieta<startEta+nCell << 133 for(int ieta=startEta; ieta<startEta+nCellInEta; ieta++){
138 for(G4int iphi=startPhi; iphi<startPhi+nCe << 134 for(int iphi=startPhi; iphi<startPhi+nCellInPhi; iphi++){
139 135
140 G4int index = this->packindex(det,z,ieta << 136 int index = this->packindex(det,z,ieta,iphi);
141 totalEnergy += themap[index]; 137 totalEnergy += themap[index];
142 #ifdef debug <<
143 ncristals+=1; 138 ncristals+=1;
>> 139 #ifdef debug
144 G4cout << "ieta - iphi - E = " << ieta < 140 G4cout << "ieta - iphi - E = " << ieta << " " << iphi << " "
145 << themap[index] << G4endl; << 141 << themap[index] << G4endl;
146 #endif 142 #endif
147 } 143 }
148 } 144 }
149 << 145
150 #ifdef debug 146 #ifdef debug
151 G4cout << "Energy in " << nCellInEta << " ce 147 G4cout << "Energy in " << nCellInEta << " cells in eta times "
152 << nCellInPhi << " cells in phi matri << 148 << nCellInPhi << " cells in phi matrix = " << totalEnergy
153 << " for " << ncristals << " crystals << 149 << " for " << ncristals << " crystals" << G4endl;
154 #endif 150 #endif
155 return totalEnergy; 151 return totalEnergy;
156 152
157 } 153 }
158 154