/[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.4 - (show annotations) (download)
Mon Mar 26 14:02:06 2007 UTC (17 years, 8 months ago) by mocchiut
Branch: MAIN
Changes since 1.3: +91 -13 lines
Bug in CaloStrip/Level0 fixed + new alignement parameters added

1 /**
2 * \file src/CaloLevel1.cpp
3 * \author Emiliano Mocchiutti
4 *
5 **/
6 #include <CaloLevel1.h>
7 //
8 ClassImp(CaloStrip);
9 ClassImp(CaloLevel1);
10
11 /**
12 * CaloStrip default constructor
13 **/
14 CaloStrip::CaloStrip() {
15 c1 = 0;
16 this->Clear();
17 UseStandardAlig();
18 };
19
20 /**
21 * CaloStrip default constructor
22 **/
23 CaloStrip::CaloStrip(CaloLevel1 *calo) {
24 c1 = calo->GetCaloLevel1();
25 this->Clear();
26 ismech = false;
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 * Connect to the DB and retrieve alignement parameters
44 **/
45 void CaloStrip::UseStandardAlig(){
46 //
47 ismech = false;
48 //
49 stringstream aligfile;
50 Int_t error = 0;
51 FILE *f = 0;
52 ifstream badfile;
53 GL_PARAM *glparam = new GL_PARAM();
54 //
55 TString host = "mysql://localhost/pamelaprod";
56 TString user = "anonymous";
57 TString psw = "";
58 const char *pamdbhost=gSystem->Getenv("PAM_DBHOST");
59 const char *pamdbuser=gSystem->Getenv("PAM_DBUSER");
60 const char *pamdbpsw=gSystem->Getenv("PAM_DBPSW");
61 if ( !pamdbhost ) pamdbhost = "";
62 if ( !pamdbuser ) pamdbuser = "";
63 if ( !pamdbpsw ) pamdbpsw = "";
64 if ( strcmp(pamdbhost,"") ) host = pamdbhost;
65 if ( strcmp(pamdbuser,"") ) user = pamdbuser;
66 if ( strcmp(pamdbpsw,"") ) psw = pamdbpsw;
67 TSQLServer *dbc = TSQLServer::Connect(host.Data(),user.Data(),psw.Data());
68 //
69 UXal = 0.;
70 UYal = 0.;
71 UZal = 0.;
72 //
73 if ( dbc ){
74 //
75 // determine where I can find calorimeter ADC to MIP conversion file
76 //
77 printf(" Querying DB for calorimeter parameters files...\n");
78 //
79 //
80 //
81 error = 0;
82 error = glparam->Query_GL_PARAM(1000,102,dbc);
83 if ( error >= 0 ){
84 //
85 aligfile.str("");
86 aligfile << glparam->PATH.Data() << "/";
87 aligfile << glparam->NAME.Data();
88 //
89 printf("\n Using parameter file: \n %s \n",aligfile.str().c_str());
90 f = fopen(aligfile.str().c_str(),"rb");
91 if ( f ){
92 //
93 fread(&UXal,sizeof(UXal),1,f);
94 fread(&UYal,sizeof(UYal),1,f);
95 fread(&UZal,sizeof(UZal),1,f);
96 //
97 fclose(f);
98 };
99 //
100 };
101 dbc->Close();
102 };
103 if ( !UXal ){
104 //
105 printf(" No able to query DB for calorimeter parameters files\n Using hard-coded parameters \n");
106 UXal = CTX;
107 UYal = CTY;
108 UZal = CTZ;
109 };
110 //
111 };
112
113
114
115
116
117
118 /**
119 * Given a strip returns its position in the PAMELA reference system
120 **/
121 void CaloStrip::Set(Int_t view, Int_t plane, Int_t strip) {
122 //
123 this->Clear();
124 //
125 if ( view < 0 || view > 1 ){
126 printf(" ERROR: 0 =< view =< 1 \n");
127 return;
128 };
129 if ( plane < 0 || plane > 21 ){
130 printf(" ERROR: 0 =< plane =< 21 \n");
131 return;
132 };
133 if ( strip < 0 || strip > 95 ){
134 printf(" ERROR: 0 =< strip =< 95 \n");
135 return;
136 };
137 //
138 Float_t lShift = 0.;
139 Float_t lPos = 0.;
140 extern struct shift shift_;
141 //
142 // Find MIPs for that strip
143 //
144 if ( c1 ) fE = c1->GetEstrip(view, plane, strip);
145 //
146 fView = view + 1;
147 fPlane = plane + 1;
148 fStrip = strip + 1;
149 if ( fPlane%2 ){
150 lShift = -0.5;
151 } else {
152 lShift = 0.5;
153 };
154 //
155 if ( fView == 2 ) lShift = - lShift;
156 //
157 shift_.shift = lShift;
158 //
159 Float_t zplane[22];
160 zplane[0] = 0.;
161 Int_t ii = 0;
162 for ( Int_t i = 2; i < 23; i++){
163 ii = i-1;
164 if ( i%2 ){
165 zplane[ii] = zplane[ii-1] - 10.09;
166 } else {
167 zplane[ii] = zplane[ii-1] - 8.09;
168 };
169 };
170 //
171 millim_(&fStrip,&lPos);
172 //
173 if ( fView == 1 ){
174 //
175 // X view
176 //
177 fX = (lPos - UXal)/10.;
178 fY = 0.;
179 fZ = (zplane[fPlane-1] - 5.81 + UZal)/10.;
180 } else {
181 //
182 // Y view
183 //
184 fX = 0;
185 fY = (lPos - UYal)/10.;
186 fZ = (zplane[fPlane-1] + UZal)/10.;
187 };
188 //
189 };
190
191 /**
192 * Given a point in the space (PAMELA ref system) returns the closest strip
193 **/
194 void CaloStrip::Set(Float_t X, Float_t Y, Float_t Z) {
195 //
196 fX = X;
197 fY = Y;
198 fZ = Z;
199 //
200 Float_t zplane[22];
201 zplane[0] = 0.;
202 Int_t ii = 0;
203 for ( Int_t i = 2; i < 23; i++){
204 ii = i-1;
205 if ( i%2 ){
206 zplane[ii] = zplane[ii-1] - 10.09;
207 } else {
208 zplane[ii] = zplane[ii-1] - 8.09;
209 };
210 };
211 //
212 Float_t dzx[22];
213 Float_t dzy[22];
214 for ( Int_t i=0; i < 22; i++){
215 dzx[i] = fabs(fZ*10. - (zplane[i] - 5.81 + UZal));
216 dzy[i] = fabs(fZ*10. - (zplane[i] + UZal));
217 };
218 //
219 Float_t minx = TMath::MinElement(22,dzx);
220 Float_t miny = TMath::MinElement(22,dzy);
221 //
222 // find view
223 //
224 if ( minx < miny ){
225 fView = 1;
226 } else {
227 fView = 2;
228 };
229 //
230 // find plane
231 //
232 Float_t pos = 0.;
233 //
234 for ( Int_t i=0; i < 22; i++){
235 if ( fView == 1 ){
236 if ( dzx[i] == minx ){
237 fPlane = i+1;
238 pos = fX*10. + UXal;
239 };
240 } else {
241 if ( dzy[i] == miny ){
242 fPlane = i+1;
243 pos = fY*10. + UYal;
244 };
245 };
246 };
247 //
248 // find strip
249 //
250 Float_t dxy[96];
251 Float_t stpos = 0.;
252 //
253 CaloStrip *ca = new CaloStrip();
254 //
255 for ( Int_t i=0; i < 96; i++){
256 ca->Set(fView-1,fPlane-1,i);
257 if ( fView == 1 ){
258 stpos = ca->GetX()*10. + UXal;
259 } else {
260 stpos = ca->GetY()*10. + UYal;
261 };
262 dxy[i] = fabs(pos - stpos);
263 };
264 //
265 delete ca;
266 //
267 Float_t mins = TMath::MinElement(96,dxy);
268 //
269 for ( Int_t i=0; i < 96; i++){
270 if ( dxy[i] == mins ) fStrip = i+1;
271 };
272 };
273
274 /**
275 * CaloLevel1 constructor
276 **/
277 CaloLevel1::CaloLevel1() {
278 //
279 estrip = TArrayI(0,NULL);
280 //
281 this->Clear();
282 //
283 };
284
285 /**
286 * Clear the CaloLevel1 object
287 **/
288 void CaloLevel1::Clear() {
289 //
290 istrip = 0;
291 estrip.Reset();
292 //
293 };
294
295 /**
296 * Returns the detected energy for the given strip once loaded the event
297 **/
298 Float_t CaloLevel1::GetEstrip(Int_t sview, Int_t splane, Int_t sstrip){
299 Int_t view = -1;
300 Int_t plane = -1;
301 Int_t strip = -1;
302 Float_t mip = 0.;
303 //
304 if ( istrip == 0 ) return(0.);
305 //
306 for (Int_t i = 0; i<istrip; i++ ){
307 //
308 mip = DecodeEstrip(i,view,plane,strip);
309 //
310 if ( view == sview && splane == plane && sstrip == strip ) return(mip);
311 //
312 // entry are ordered by strip, plane and view number. Go out if you pass the input strip
313 //
314 if ( view == sview && plane == splane && strip > sstrip ) return(0.);
315 if ( view == sview && plane > splane ) return(0.);
316 if ( view > sview ) return(0.);
317 //
318 };
319 return(0.);
320 };
321
322 /**
323 * Given estrip entry returns energy plus view, plane and strip numbers
324 **/
325 Float_t CaloLevel1::DecodeEstrip(Int_t entry, Int_t &view, Int_t &plane, Int_t &strip){
326 //
327 if ( entry>istrip ) return(0.);
328 //
329 // printf(" num lim %f \n",std::numeric_limits<Float_t>::max());
330 // printf(" estrip.At(%i) = %i \n",entry,estrip.At(entry));
331 //
332 Int_t eval = 0;
333 if ( estrip.At(entry) > 0. ){
334 view = 0;
335 eval = estrip.At(entry);
336 } else {
337 view = 1;
338 eval = -estrip.At(entry);
339 };
340 //
341 Int_t fbi = 0;
342 fbi = (Int_t)truncf((Float_t)(eval/1000000000));
343 //
344 Int_t plom = 0;
345 plom = (Int_t)truncf((Float_t)((eval-fbi*1000000000)/10000000));
346 //
347 Float_t tim = 100000.;
348 plane = plom;
349 if ( fbi == 1 ) tim = 10000.;
350 if ( plom > 21 ){
351 plane = plom - 22;
352 if ( fbi == 1 ){
353 tim = 1000.;
354 } else {
355 tim = 100.;
356 };
357 };
358 if ( plom > 43 ){
359 plane = plom - 44;
360 tim = 10.;
361 };
362 if ( plom > 65 ){
363 plane = plom - 66;
364 tim = 1.;
365 };
366 //
367 strip = (Int_t)truncf((Float_t)((eval - fbi*1000000000 -plom*10000000)/100000));
368 //
369 Float_t mip = ((Float_t)(eval - fbi*1000000000 -plom*10000000 -strip*100000))/tim;
370 //
371 if ( mip > 0. && mip < 99999. ) return(mip);
372 //
373 printf(" WARNING: problems decoding value %i at entry %i \n",estrip.At(entry),entry);
374 //
375 view = -1;
376 plane = -1;
377 strip = -1;
378 return(0.);
379 }
380
381 /*
382 * Returns energy released on plane nplane (where 0<= nplane <= 43, 0 = 1Y, 1 = 1X, 2 = 2Y, 3 = 2X, etc. etc.).
383 */
384 Float_t CaloLevel1::qtotpl(Int_t nplane){
385 //
386 Int_t sview = 1;
387 if ( nplane%2 ) sview = 0;
388 //
389 Int_t splane = nplane-(sview+1)/2;
390 //
391 Float_t totmip = qtotpl(sview,splane);
392 //
393 return(totmip);
394 //
395 };
396
397 /*
398 * Returns energy released on view "view" (0 = X, 1 = Y) and plane "plane" ( 0 <= plane <= 21 ).
399 */
400 Float_t CaloLevel1::qtotpl(Int_t sview, Int_t splane){
401 //
402 Int_t view = -1;
403 Int_t plane = -1;
404 Int_t strip = -1;
405 //
406 Float_t mip = 0.;
407 Float_t totmip = 0.;
408 //
409 if ( istrip == 0 ) return(0.);
410 //
411 for (Int_t i = 0; i<istrip; i++ ){
412 //
413 mip = DecodeEstrip(i,view,plane,strip);
414 //
415 if ( view == sview && splane == plane ) totmip += mip;
416 //
417 // entry are ordered by strip, plane and view number. Go out if you pass the input strip
418 //
419 if ( view == sview && plane > splane ) return(totmip);
420 if ( view > sview ) return(totmip);
421 //
422 };
423 //
424 return(totmip);
425 //
426 };

  ViewVC Help
Powered by ViewVC 1.1.23