--- calo/flight/CaloEnergy/src/CaloEnergy.cpp 2009/08/11 13:28:27 1.5 +++ calo/flight/CaloEnergy/src/CaloEnergy.cpp 2009/08/17 15:00:53 1.11 @@ -46,6 +46,7 @@ // cp->ForceCaloFit(); // cp->SetDebug(true); // cp->Process(); + if ( clong ) clong->SetCaloLevel2Pointer(cp->GetLevel2Pointer()); } @@ -189,8 +190,10 @@ // Use CaloStrip to determine once the position of border strips for each section // // - fM = 2. + 0.096; // real position from cbar - fM1 = 2. - 0.122 - 0.096; // due to calculation of xe1 etc. + // fM = 2. + 0.096; // real position from cbar BUG the 0.096 is already taken into account in the border calculation made by Giovanna + fM = 2. ; // real position from cbar + // fM1 = 2. - 0.122 - 0.096; // due to calculation of xe1 etc. BUG! this way we count from the silicon border not from the silicon sensitive area + fM1 = 2. - 0.122 - 0.096 + 0.096; // due to calculation of xe1 etc. if ( fM1 < 0. ) fM1 = 0.; // CaloStrip *cs = new CaloStrip(fSimu); @@ -215,6 +218,8 @@ cs->Set(1,0,95); xe6= cs->GetY(); + // printf(" xe1 %f xe2 %f xe3 %f xe4 %f xe5 %f xe6 %f \n",xe1,xe2,xe3,xe4,xe5,xe6); + // view x plane 0 strip 0 cs->Set(0,0,0); yo1= cs->GetX(); @@ -235,7 +240,6 @@ cs->Set(0,0,95); yo6= cs->GetX(); - // view y plane 1 strip 0 cs->Set(1,1,0); xo1= cs->GetY(); @@ -276,6 +280,8 @@ cs->Set(0,1,95); ye6= cs->GetX(); + //printf(" ye1 %f ye2 %f ye3 %f ye4 %f ye5 %f ye6 %f \n",ye1,ye2,ye3,ye4,ye5,ye6); + // for (Int_t p = 0; p<22; p ++){ for (Int_t v = 0; v<2; v++ ){ @@ -314,7 +320,7 @@ yomax_en= 0.; yemax_en= 0.; // - memset(enstrip,0,2*22*96*(sizeof(Float_t))); + memset(enstrip,0,2*22*96*(sizeof(Float_t))); en = 0.; view = 0; plane = 0; @@ -377,6 +383,10 @@ printf(" energyxo:.............. %f \n",energyxo); printf(" energyye:.............. %f \n",energyye); printf(" energyyo:.............. %f \n",energyyo); + printf(" fXEen_maxplane:........ %f \n",fXEen_maxplane); + printf(" fXOen_maxplane:........ %f \n",fXOen_maxplane); + printf(" fYEen_maxplane:........ %f \n",fYEen_maxplane); + printf(" fYOen_maxplane:........ %f \n",fYOen_maxplane); printf(" x0max :.............. %f \n",x0max); printf(" debug :.............. %i \n",debug); @@ -453,13 +463,15 @@ enstrip[view][plane][strip]=en; }; // + if ( debug && ((fM1+0.122-0.244*(Float_t)fRad) < 0.) ) printf("Error: (fM1+0.122-0.244*(Float_t)fRad) < 0. fM1 %f fRad %i %f \n",fM1,fRad,(fM1+0.122-0.244*(Float_t)fRad)); + // // sum energy plane by plane for each sections // for (Int_t i=0;i<11;i++){ for(strip=0; strip<96; strip++) { if ( fRad < 0 ){ // - // run over all the strips of the plane + // run over all the strips of the plane // en_xep[i] += enstrip[1][2*i][strip]; en_yop[i] += enstrip[0][2*i][strip]; @@ -469,10 +481,10 @@ // // use only the strips inside a cylinder of given radius fRad // - if ( strip >= cl2->cibar[2*i][1]-1-fRad && strip <= cl2->cibar[2*i][1]-1+fRad ) en_xep[i] += enstrip[1][2*i][strip]; - if ( strip >= cl2->cibar[2*i][0]-1-fRad && strip <= cl2->cibar[2*i][0]-1+fRad ) en_yop[i] += enstrip[0][2*i][strip]; - if ( strip >= cl2->cibar[(2*i)+1][1]-1-fRad && strip <= cl2->cibar[(2*i)+1][1]-1+fRad ) en_xop[i] += enstrip[1][(2*i)+1][strip]; - if ( strip >= cl2->cibar[(2*i)+1][0]-1-fRad && strip <= cl2->cibar[(2*i)+1][0]-1+fRad ) en_yep[i] += enstrip[0][(2*i)+1][strip]; + if ( cl2->cibar[2*i][1] >= 1 && cl2->cibar[2*i][1] <= 96 && (strip >= (cl2->cibar[2*i][1]-1-fRad)) && (strip <= (cl2->cibar[2*i][1]-1+fRad)) ) en_xep[i] += enstrip[1][2*i][strip]; + if ( cl2->cibar[2*i][0] >= 1 && cl2->cibar[2*i][0] <= 96 && (strip >= (cl2->cibar[2*i][0]-1-fRad)) && (strip <= (cl2->cibar[2*i][0]-1+fRad)) ) en_yop[i] += enstrip[0][2*i][strip]; + if ( cl2->cibar[(2*i)+1][1] >= 1 && cl2->cibar[(2*i)+1][1] <= 96 && (strip >= (cl2->cibar[(2*i)+1][1]-1-fRad)) && (strip <= (cl2->cibar[(2*i)+1][1]-1+fRad)) ) en_xop[i] += enstrip[1][(2*i)+1][strip]; + if ( cl2->cibar[(2*i)+1][0] >= 1 && cl2->cibar[(2*i)+1][0] <= 96 && (strip >= (cl2->cibar[(2*i)+1][0]-1-fRad)) && (strip <= (cl2->cibar[(2*i)+1][0]-1+fRad)) ) en_yep[i] += enstrip[0][(2*i)+1][strip]; }; }; energyxe += en_xep[i]; @@ -538,38 +550,47 @@ Int_t wpl = (Int_t)roundf(x0max/0.76); Bool_t isY = false; if ( ((x0max/0.76)-(Float_t)wpl) > 0. ) isY = true; + xomax_en = 0.; + yemax_en = 0.; + xemax_en = 0.; + yomax_en = 0.; + // if ( !(wpl%2) ){ // 0, 2, 4, ... if ( isY ){ - yomax_en = 1000.; - xemax_en = 500.; + if ( section.Contains("YO") ) yomax_en = 1000.; + if ( section.Contains("XE") ) xemax_en = 500.; fMax_planeyo=wpl/2; fMax_planexe=wpl/2; + if ( section.Contains("XO") ) xomax_en = 10.; + if ( section.Contains("YE") ) yemax_en = 5.; } else { - yomax_en = 500.; - xemax_en = 1000.; + if ( section.Contains("YO") ) yomax_en = 500.; + if ( section.Contains("XE") ) xemax_en = 1000.; fMax_planeyo=wpl/2; fMax_planexe=wpl/2; + if ( section.Contains("XO") ) xomax_en = 5.; + if ( section.Contains("YE") ) yemax_en = 10.; }; - xomax_en = 0.; - yemax_en = 0.; } else { // 1, 3, 5, ... if ( isY ){ - yemax_en = 1000.; - xomax_en = 500.; + if ( section.Contains("YE") ) yemax_en = 1000.; + if ( section.Contains("XO") ) xomax_en = 500.; fMax_planeye=(wpl-1)/2; fMax_planexo=(wpl-1)/2; + if ( section.Contains("XE") ) xemax_en = 10.; + if ( section.Contains("YO") ) yomax_en = 5.; } else { - yemax_en = 500.; - xomax_en = 1000.; + if ( section.Contains("YE") ) yemax_en = 500.; + if ( section.Contains("XO") ) xomax_en = 1000.; fMax_planeye=(wpl-1)/2; fMax_planexo=(wpl-1)/2; + if ( section.Contains("XE") ) xemax_en = 5.; + if ( section.Contains("YO") ) yomax_en = 10.; }; - xemax_en = 0.; - yomax_en = 0.; }; - if ( debug ) printf(" x0max %f wpl %i isY %i yomax_en %f xemax_en %f yemax_en %f xomax_en %f fMaxplane %i %i %i %i\n",x0max,wpl,isY,yomax_en,xemax_en,yemax_en,xomax_en,fMax_planeyo,fMax_planexe,fMax_planeye,fMax_planexo); + if ( debug ) printf(" x0max %f x0max/0.76 %f wpl %i isY %i yomax_en %f xemax_en %f yemax_en %f xomax_en %f fMaxplane %i %i %i %i\n",x0max,(x0max/0.76),wpl,isY,yomax_en,xemax_en,yemax_en,xomax_en,fMax_planeyo,fMax_planexe,fMax_planeye,fMax_planexo); }; // Int_t nPl = fPl; @@ -641,8 +662,6 @@ // // for each plane of the calorimeter find the position of the track in the direction along the strip (where we do not have a measurement from the selected plane) by looking at the plane above/below of the other view and extrapolating the trajectory to the given plane // - Float_t track_coordx[22][2]; - Float_t track_coordy[22][2]; // Float_t tgx_cl2; Float_t tgy_cl2; @@ -652,8 +671,10 @@ for (Int_t p=0; p<22; p++){ track_coordy[p][1] = cl2->cbar[p][1]; track_coordx[p][1] = cl2->cbar[p][0] - fabs(trk_z[p][1]-trk_z[p][0])*tgx_cl2; + // track_coordx[p][1] = cl2->cbar[p][0] + fabs(trk_z[p][1]-trk_z[p][0])*tgx_cl2; track_coordx[p][0] = cl2->cbar[p][0]; track_coordy[p][0] = cl2->cbar[p][1] - fabs(trk_z[p][1]-trk_z[p][0])*tgy_cl2; + // track_coordy[p][0] = cl2->cbar[p][1] + fabs(trk_z[p][1]-trk_z[p][0])*tgy_cl2; if ( debug ) printf(" p %i track_coordy[p][1] %f track_coordx[p][1] %f track_coordx[p][0] %f track_coordy[p][0] %f \n",p,track_coordy[p][1],track_coordx[p][1],track_coordx[p][0],track_coordy[p][0]); }; // @@ -712,7 +733,10 @@ // event is contained (or partially contained) hence we can integrate energy up to the maximum and calculate the energy as measured by this section // if ( fXosel ){ - for (Int_t iplm=0;iplm<=TMath::Min(10,(Int_t)(fMax_planexo+nPl)) ;iplm++) fXOen_maxplane += en_xop[iplm]; + for (Int_t iplm=0;iplm<=TMath::Min(10,(Int_t)(fMax_planexo+nPl)) ;iplm++){ + fXOen_maxplane += en_xop[iplm]; + if ( debug ) printf(" XO iplm %i fXOen_maxplane %f en_xop[iplm] %f\n",iplm,fXOen_maxplane,en_xop[iplm]); + }; fEnergyxo = fXOen_maxplane/fConv_rxo; }; }; @@ -761,7 +785,10 @@ fXesel = true; }; if ( fXesel ){ - for (Int_t iplm=0;iplm<=TMath::Min(10,(Int_t)(fMax_planexe+nPl)) ;iplm++) fXEen_maxplane += en_xep[iplm]; + for (Int_t iplm=0;iplm<=TMath::Min(10,(Int_t)(fMax_planexe+nPl)) ;iplm++){ + fXEen_maxplane += en_xep[iplm]; + if ( debug ) printf(" XE iplm %i fXOen_maxplane %f en_xop[iplm] %f\n",iplm,fXEen_maxplane,en_xep[iplm]); + }; fEnergyxe = fXEen_maxplane/fConv_rxe; }; }; @@ -971,7 +998,12 @@ } else { clong->Fit(); }; - fXOen_maxplane = clong->Get_E0(); + if ( clong->GetLowerLimit() != 0. || clong->GetUpperLimit() != 0. ){ + fXOen_maxplane = clong->Get_defE0(); + } else { + fXOen_maxplane = clong->Get_E0(); + }; + fMax_plane = clong->Get_tmax(); fYOen_maxplane = 0.; fYEen_maxplane = 0.; fXEen_maxplane = 0.;