/[PAMELA software]/DarthVader/CalorimeterLevel2/src/CaloLevel1.cpp
ViewVC logotype

Contents of /DarthVader/CalorimeterLevel2/src/CaloLevel1.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.3 - (show annotations) (download)
Thu Jan 11 16:08:02 2007 UTC (17 years, 10 months ago) by mocchiut
Branch: MAIN
CVS Tags: v3r00
Changes since 1.2: +46 -0 lines
Calorimeter: added methods to get released energy on any plane

1 /**
2 * \file src/CaloLevel1.cpp
3 * \author Emiliano Mocchiutti
4 *
5 **/
6 #include <TObject.h>
7 #include <TMath.h>
8 #include <CaloLevel1.h>
9 //
10 ClassImp(CaloStrip);
11 ClassImp(CaloLevel1);
12
13 /**
14 * CaloStrip default constructor
15 **/
16 CaloStrip::CaloStrip() {
17 c1 = 0;
18 this->Clear();
19 };
20
21 /**
22 * CaloStrip default constructor
23 **/
24 CaloStrip::CaloStrip(CaloLevel1 *calo) {
25 c1 = calo->GetCaloLevel1();
26 this->Clear();
27 };
28
29 /**
30 * Clear variables
31 **/
32 void CaloStrip::Clear() {
33 fE = 0.;
34 fX = 0.;
35 fY = 0.;
36 fZ = 0.;
37 fView = 0;
38 fPlane = 0;
39 fStrip = 0;
40 };
41
42 /**
43 * Given a strip returns its position in the PAMELA reference system
44 **/
45 void CaloStrip::Set(Int_t view, Int_t plane, Int_t strip) {
46 //
47 this->Clear();
48 //
49 if ( view < 0 || view > 1 ){
50 printf(" ERROR: 0 =< view =< 1 \n");
51 return;
52 };
53 if ( plane < 0 || plane > 21 ){
54 printf(" ERROR: 0 =< plane =< 21 \n");
55 return;
56 };
57 if ( strip < 0 || strip > 95 ){
58 printf(" ERROR: 0 =< strip =< 95 \n");
59 return;
60 };
61 //
62 Float_t lShift = 0.;
63 Float_t lPos = 0.;
64 extern struct shift shift_;
65 //
66 // Find MIPs for that strip
67 //
68 if ( c1 ) fE = c1->GetEstrip(view, plane, strip);
69 //
70 fView = view + 1;
71 fPlane = plane + 1;
72 fStrip = strip + 1;
73 if ( fPlane%2 ){
74 lShift = -0.5;
75 } else {
76 lShift = 0.5;
77 };
78 //
79 shift_.shift = lShift;
80 //
81 Float_t zplane[22];
82 zplane[0] = 0.;
83 Int_t ii = 0;
84 for ( Int_t i = 2; i < 23; i++){
85 ii = i-1;
86 if ( i%2 ){
87 zplane[ii] = zplane[ii-1] - 10.09;
88 } else {
89 zplane[ii] = zplane[ii-1] - 8.09;
90 };
91 };
92 //
93 millim_(&fStrip,&lPos);
94 //
95 if ( fView == 1 ){
96 //
97 // X view
98 //
99 fX = (lPos - CTX)/10.;
100 fY = 0.;
101 fZ = (zplane[fPlane-1] - 5.81 + CTZ)/10.;
102 } else {
103 //
104 // Y view
105 //
106 fX = 0;
107 fY = (lPos - CTY)/10.;
108 fZ = (zplane[fPlane-1] + CTZ)/10.;
109 };
110 //
111 };
112
113 /**
114 * Given a point in the space (PAMELA ref system) returns the closest strip
115 **/
116 void CaloStrip::Set(Float_t X, Float_t Y, Float_t Z) {
117 //
118 fX = X;
119 fY = Y;
120 fZ = Z;
121 //
122 Float_t zplane[22];
123 zplane[0] = 0.;
124 Int_t ii = 0;
125 for ( Int_t i = 2; i < 23; i++){
126 ii = i-1;
127 if ( i%2 ){
128 zplane[ii] = zplane[ii-1] - 10.09;
129 } else {
130 zplane[ii] = zplane[ii-1] - 8.09;
131 };
132 };
133 //
134 Float_t dzx[22];
135 Float_t dzy[22];
136 for ( Int_t i=0; i < 22; i++){
137 dzx[i] = fabs(fZ*10. - (zplane[i] - 5.81 + CTZ));
138 dzy[i] = fabs(fZ*10. - (zplane[i] + CTZ));
139 };
140 //
141 Float_t minx = TMath::MinElement(22,dzx);
142 Float_t miny = TMath::MinElement(22,dzy);
143 //
144 // find view
145 //
146 if ( minx < miny ){
147 fView = 1;
148 } else {
149 fView = 2;
150 };
151 //
152 // find plane
153 //
154 Float_t pos = 0.;
155 //
156 for ( Int_t i=0; i < 22; i++){
157 if ( fView == 1 ){
158 if ( dzx[i] == minx ){
159 fPlane = i+1;
160 pos = fX*10. + CTX;
161 };
162 } else {
163 if ( dzy[i] == miny ){
164 fPlane = i+1;
165 pos = fY*10. + CTY;
166 };
167 };
168 };
169 //
170 // find strip
171 //
172 Float_t dxy[96];
173 Float_t stpos = 0.;
174 //
175 CaloStrip *ca = new CaloStrip();
176 //
177 for ( Int_t i=0; i < 96; i++){
178 ca->Set(fView-1,fPlane-1,i);
179 if ( fView == 1 ){
180 stpos = ca->GetX()*10. + CTX;
181 } else {
182 stpos = ca->GetY()*10. + CTY;
183 };
184 dxy[i] = fabs(pos - stpos);
185 };
186 //
187 delete ca;
188 //
189 Float_t mins = TMath::MinElement(96,dxy);
190 //
191 for ( Int_t i=0; i < 96; i++){
192 if ( dxy[i] == mins ) fStrip = i+1;
193 };
194 };
195
196 /**
197 * CaloLevel1 constructor
198 **/
199 CaloLevel1::CaloLevel1() {
200 //
201 estrip = TArrayI(0,NULL);
202 //
203 this->Clear();
204 //
205 };
206
207 /**
208 * Clear the CaloLevel1 object
209 **/
210 void CaloLevel1::Clear() {
211 //
212 istrip = 0;
213 estrip.Reset();
214 //
215 };
216
217 /**
218 * Returns the detected energy for the given strip once loaded the event
219 **/
220 Float_t CaloLevel1::GetEstrip(Int_t sview, Int_t splane, Int_t sstrip){
221 Int_t view = -1;
222 Int_t plane = -1;
223 Int_t strip = -1;
224 Float_t mip = 0.;
225 //
226 if ( istrip == 0 ) return(0.);
227 //
228 for (Int_t i = 0; i<istrip; i++ ){
229 //
230 mip = DecodeEstrip(i,view,plane,strip);
231 //
232 if ( view == sview && splane == plane && sstrip == strip ) return(mip);
233 //
234 // entry are ordered by strip, plane and view number. Go out if you pass the input strip
235 //
236 if ( view == sview && plane == splane && strip > sstrip ) return(0.);
237 if ( view == sview && plane > splane ) return(0.);
238 if ( view > sview ) return(0.);
239 //
240 };
241 return(0.);
242 };
243
244 /**
245 * Given estrip entry returns energy plus view, plane and strip numbers
246 **/
247 Float_t CaloLevel1::DecodeEstrip(Int_t entry, Int_t &view, Int_t &plane, Int_t &strip){
248 //
249 if ( entry>istrip ) return(0.);
250 //
251 // printf(" num lim %f \n",std::numeric_limits<Float_t>::max());
252 // printf(" estrip.At(%i) = %i \n",entry,estrip.At(entry));
253 //
254 Int_t eval = 0;
255 if ( estrip.At(entry) > 0. ){
256 view = 0;
257 eval = estrip.At(entry);
258 } else {
259 view = 1;
260 eval = -estrip.At(entry);
261 };
262 //
263 Int_t fbi = 0;
264 fbi = (Int_t)truncf((Float_t)(eval/1000000000));
265 //
266 Int_t plom = 0;
267 plom = (Int_t)truncf((Float_t)((eval-fbi*1000000000)/10000000));
268 //
269 Float_t tim = 100000.;
270 plane = plom;
271 if ( fbi == 1 ) tim = 10000.;
272 if ( plom > 21 ){
273 plane = plom - 22;
274 if ( fbi == 1 ){
275 tim = 1000.;
276 } else {
277 tim = 100.;
278 };
279 };
280 if ( plom > 43 ){
281 plane = plom - 44;
282 tim = 10.;
283 };
284 if ( plom > 65 ){
285 plane = plom - 66;
286 tim = 1.;
287 };
288 //
289 strip = (Int_t)truncf((Float_t)((eval - fbi*1000000000 -plom*10000000)/100000));
290 //
291 Float_t mip = ((Float_t)(eval - fbi*1000000000 -plom*10000000 -strip*100000))/tim;
292 //
293 if ( mip > 0. && mip < 99999. ) return(mip);
294 //
295 printf(" WARNING: problems decoding value %i at entry %i \n",estrip.At(entry),entry);
296 //
297 view = -1;
298 plane = -1;
299 strip = -1;
300 return(0.);
301 }
302
303 /*
304 * Returns energy released on plane nplane (where 0<= nplane <= 43, 0 = 1Y, 1 = 1X, 2 = 2Y, 3 = 2X, etc. etc.).
305 */
306 Float_t CaloLevel1::qtotpl(Int_t nplane){
307 //
308 Int_t sview = 1;
309 if ( nplane%2 ) sview = 0;
310 //
311 Int_t splane = nplane-(sview+1)/2;
312 //
313 Float_t totmip = qtotpl(sview,splane);
314 //
315 return(totmip);
316 //
317 };
318
319 /*
320 * Returns energy released on view "view" (0 = X, 1 = Y) and plane "plane" ( 0 <= plane <= 21 ).
321 */
322 Float_t CaloLevel1::qtotpl(Int_t sview, Int_t splane){
323 //
324 Int_t view = -1;
325 Int_t plane = -1;
326 Int_t strip = -1;
327 //
328 Float_t mip = 0.;
329 Float_t totmip = 0.;
330 //
331 if ( istrip == 0 ) return(0.);
332 //
333 for (Int_t i = 0; i<istrip; i++ ){
334 //
335 mip = DecodeEstrip(i,view,plane,strip);
336 //
337 if ( view == sview && splane == plane ) totmip += mip;
338 //
339 // entry are ordered by strip, plane and view number. Go out if you pass the input strip
340 //
341 if ( view == sview && plane > splane ) return(totmip);
342 if ( view > sview ) return(totmip);
343 //
344 };
345 //
346 return(totmip);
347 //
348 };

  ViewVC Help
Powered by ViewVC 1.1.23