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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.8 - (hide annotations) (download)
Fri Mar 30 14:37:44 2007 UTC (17 years, 8 months ago) by mocchiut
Branch: MAIN
CVS Tags: v3r04, v3r05, v3r03, v4r00, v3r06
Changes since 1.7: +64 -52 lines
CaloStrip static alignement variables (no DB connection at any time)

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

  ViewVC Help
Powered by ViewVC 1.1.23