From: pavlinov <pavlinov@f7af4fe6-9843-0410-8265-dc069ae4e863>
Date: Tue, 11 Sep 2007 19:43:55 +0000 (+0000)
Subject: minor corrections for calibration staff
X-Git-Url: http://git.uio.no/git/?p=u%2Fmrichter%2FAliRoot.git;a=commitdiff_plain;h=1ae500a2aac37a680262f44413990e407a73f47f

minor corrections for calibration staff
---

diff --git a/EMCAL/AliEMCALGeometry.cxx b/EMCAL/AliEMCALGeometry.cxx
index 32ed5c939ba..f77109649d4 100644
--- a/EMCAL/AliEMCALGeometry.cxx
+++ b/EMCAL/AliEMCALGeometry.cxx
@@ -47,6 +47,7 @@
 #include <TMatrixD.h>
 #include <TObjArray.h>
 #include <TObjString.h>
+#include <TVector2.h>
 #include <TVector3.h>
 
 // -- ALICE Headers.
@@ -334,7 +335,7 @@ void AliEMCALGeometry::Init(void){
       fECScintThick  = fECPbRadThickness = 0.5;
     }
     if(fGeoName.Contains("WSUC")){ // 18-may-05 - about common structure
-      fShellThickness = 30.; // should be change 
+      fShellThickness = 30.;       // should be change 
       fNPhi = fNZ = 4; 
     }
 
@@ -439,6 +440,9 @@ void AliEMCALGeometry::Init(void){
 					   xpos,ypos,-zpos, geoRot1);
     } // for
   }
+
+  if(fGeoName.Contains("WSUC")) fNumberOfSuperModules = 1; // Jul 12, 2007
+
   fgInit = kTRUE; 
   AliInfo(" is ended");  
 }
@@ -934,6 +938,110 @@ Bool_t AliEMCALGeometry::RelPosCellInSModule(Int_t absId, TVector3 &vloc) const
   // Alice numbering scheme - Jun 03, 2006
 }
 
+Bool_t AliEMCALGeometry::RelPosCellInSModule(Int_t absId, Double_t distEff, Double_t &xr, Double_t &yr, Double_t &zr) const
+{
+  // Jul 30, 2007 - taking into account position of shower max
+  // Look to see what the relative
+  // position inside a given cell is
+  // for a recpoint.
+  // In:
+  // absId   - cell is as in Geant,     0<= absId   < fNCells;
+  // e       - cluster energy
+  // OUT:
+  // xr,yr,zr - x,y,z coordinates of cell with absId inside SM 
+
+  // Shift index taking into account the difference between standard SM 
+  // and SM of half size in phi direction
+  const  Int_t phiIndexShift = fCentersOfCellsPhiDir.GetSize()/4; // Nov 22, 2006; was 6 for cas 2X2
+  static Int_t nSupMod, nModule, nIphi, nIeta, iphi, ieta;
+  static Int_t iphim, ietam;
+  static AliEMCALShishKebabTrd1Module *mod = 0;
+  static TVector2 v;
+  if(!CheckAbsCellId(absId)) return kFALSE;
+
+  GetCellIndex(absId, nSupMod, nModule, nIphi, nIeta);
+  GetModulePhiEtaIndexInSModule(nSupMod, nModule, iphim, ietam);
+  GetCellPhiEtaIndexInSModule(nSupMod,nModule,nIphi,nIeta, iphi, ieta); 
+ 
+  mod = GetShishKebabModule(ietam);
+  mod->GetPositionAtCenterCellLine(nIeta, distEff, v); 
+  xr = v.Y() - fParSM[0];
+  zr = v.X() - fParSM[2];
+
+  if(nSupMod<10) {
+    yr = fCentersOfCellsPhiDir.At(iphi);
+  } else {
+    yr = fCentersOfCellsPhiDir.At(iphi + phiIndexShift);
+  }
+  AliDebug(1,Form("absId %i nSupMod %i iphi %i ieta %i xr %f yr %f zr %f ",absId,nSupMod,iphi,ieta,xr,yr,zr));
+
+  return kTRUE;
+}
+
+Bool_t AliEMCALGeometry::RelPosCellInSModule(Int_t absId, Int_t maxAbsId, Double_t distEff, Double_t &xr, Double_t &yr, Double_t &zr) const
+{
+  // Jul 31, 2007 - taking into account position of shower max and apply coor2.
+  // Look to see what the relative
+  // position inside a given cell is
+  // for a recpoint.
+  // In:
+  // absId     - cell is as in Geant,     0<= absId   < fNCells;
+  // maxAbsId  - abs id of cell with highest energy
+  // e         - cluster energy
+  // OUT:
+  // xr,yr,zr - x,y,z coordinates of cell with absId inside SM 
+
+  // Shift index taking into account the difference between standard SM 
+  // and SM of half size in phi direction
+  const  Int_t phiIndexShift = fCentersOfCellsPhiDir.GetSize()/4; // Nov 22, 2006; was 6 for cas 2X2
+  static Int_t nSupMod, nModule, nIphi, nIeta, iphi, ieta;
+  static Int_t iphim, ietam;
+  static AliEMCALShishKebabTrd1Module *mod = 0;
+  static TVector2 v;
+
+  static Int_t nSupModM, nModuleM, nIphiM, nIetaM, iphiM, ietaM;
+  static Int_t iphimM, ietamM, maxAbsIdCopy=-1;
+  static AliEMCALShishKebabTrd1Module *modM = 0;
+  static Double_t distCorr;
+
+  if(!CheckAbsCellId(absId)) return kFALSE;
+
+  GetCellIndex(absId, nSupMod, nModule, nIphi, nIeta);
+  GetModulePhiEtaIndexInSModule(nSupMod, nModule, iphim, ietam);
+  GetCellPhiEtaIndexInSModule(nSupMod,nModule,nIphi,nIeta, iphi, ieta); 
+  mod = GetShishKebabModule(ietam);
+
+  if(absId != maxAbsId) {
+    distCorr = 0.;
+    if(maxAbsIdCopy != maxAbsId) {
+      GetCellIndex(maxAbsId, nSupModM, nModuleM, nIphiM, nIetaM);
+      GetModulePhiEtaIndexInSModule(nSupModM, nModuleM, iphimM, ietamM);
+      GetCellPhiEtaIndexInSModule(nSupModM,nModuleM,nIphiM,nIetaM, iphiM, ietaM); 
+      modM = GetShishKebabModule(ietamM); // do I need this ?
+      maxAbsIdCopy = maxAbsId;
+    }
+    if(ietamM !=0) {
+      distCorr = GetEtaModuleSize()*(ietam-ietamM)/TMath::Tan(modM->GetTheta()); // Stay here
+      //printf(" distCorr %f | dist %f | ietam %i -> etamM %i\n", distCorr, dist, ietam, ietamM);  
+    }
+    // distEff += distCorr;
+  }
+  // Bad resolution in this case, strong bias vs phi
+  // distEff = 0.0; 
+  mod->GetPositionAtCenterCellLine(nIeta, distEff, v); // Stay here
+  xr = v.Y() - fParSM[0];
+  zr = v.X() - fParSM[2];
+
+  if(nSupMod<10) {
+    yr = fCentersOfCellsPhiDir.At(iphi);
+  } else {
+    yr = fCentersOfCellsPhiDir.At(iphi + phiIndexShift);
+  }
+  AliDebug(1,Form("absId %i nSupMod %i iphi %i ieta %i xr %f yr %f zr %f ",absId,nSupMod,iphi,ieta,xr,yr,zr));
+
+  return kTRUE;
+}
+
 void AliEMCALGeometry::CreateListOfTrd1Modules()
 {
   // Generate the list of Trd1 modules
@@ -1268,7 +1376,7 @@ Bool_t AliEMCALGeometry::GetAbsCellIdFromEtaPhi(Double_t eta, Double_t phi, Int_
   return kFALSE;
 }
 
-AliEMCALShishKebabTrd1Module* AliEMCALGeometry::GetShishKebabModule(Int_t neta)
+AliEMCALShishKebabTrd1Module* AliEMCALGeometry::GetShishKebabModule(Int_t neta) const
 {
   //This method was too long to be
   //included in the header file - the
diff --git a/EMCAL/AliEMCALGeometry.h b/EMCAL/AliEMCALGeometry.h
index 83085db08e9..529ee219824 100644
--- a/EMCAL/AliEMCALGeometry.h
+++ b/EMCAL/AliEMCALGeometry.h
@@ -142,7 +142,7 @@ public:
   // TRD1 staff
   void    CreateListOfTrd1Modules();
   TList  *GetShishKebabTrd1Modules() const {return fShishKebabTrd1Modules;}
-  AliEMCALShishKebabTrd1Module *GetShishKebabModule(Int_t neta);
+  AliEMCALShishKebabTrd1Module *GetShishKebabModule(Int_t neta) const;
 
   void     GetTransformationForSM();
   Float_t *GetSuperModulesPars() {return fParSM;}
@@ -190,6 +190,10 @@ public:
   Bool_t   RelPosCellInSModule(Int_t absId, Double_t &xr, Double_t &yr, Double_t &zr) const;
   Bool_t   RelPosCellInSModule(Int_t absId, Double_t loc[3]) const;
   Bool_t   RelPosCellInSModule(Int_t absId, TVector3 &vloc) const;
+  //  Methods for AliEMCALRecPoint with taking into account energy of rec.point - Jul 30. 2007
+  Bool_t RelPosCellInSModule(Int_t absId,Double_t distEff,Double_t &xr,Double_t &yr,Double_t & zr) const; 
+  Bool_t RelPosCellInSModule(Int_t absId,Int_t maxAbsId,Double_t distEff,Double_t &xr,Double_t &yr,Double_t &zr) const;
+
   // ---
   Float_t AngleFromEta(Float_t eta) const { // returns theta in radians for a given pseudorapidity
     return 2.0*TMath::ATan(TMath::Exp(-eta));
diff --git a/EMCAL/AliEMCALRecPoint.cxx b/EMCAL/AliEMCALRecPoint.cxx
index 3503dfccd66..c0b8207f94f 100644
--- a/EMCAL/AliEMCALRecPoint.cxx
+++ b/EMCAL/AliEMCALRecPoint.cxx
@@ -61,7 +61,8 @@ AliEMCALRecPoint::AliEMCALRecPoint()
     fMaxParent(0),
     fParentsList(0),
     fDEParentsList(0),
-    fSuperModuleNumber(0)
+    fSuperModuleNumber(0),
+    fDigitIndMax(-1)
 {
   // ctor
   AliRunLoader *rl = AliRunLoader::GetRunLoader();
@@ -90,7 +91,8 @@ AliEMCALRecPoint::AliEMCALRecPoint(const char * opt)
     fMaxParent(1000),
     fParentsList(0),
     fDEParentsList(0),
-    fSuperModuleNumber(0)
+    fSuperModuleNumber(0),
+    fDigitIndMax(-1)
 {
   // ctor
   // Increase fMaxTrack for EMCAL.
@@ -134,7 +136,8 @@ AliEMCALRecPoint::AliEMCALRecPoint(const AliEMCALRecPoint & rp)
     fMaxParent(rp.fMaxParent),
     fParentsList(0),
     fDEParentsList(0),
-    fSuperModuleNumber(rp.fSuperModuleNumber)
+    fSuperModuleNumber(rp.fSuperModuleNumber),
+    fDigitIndMax(rp.fDigitIndMax)
 {
   //copy ctor
   fLambda[0] = rp.fLambda[0];
@@ -470,15 +473,25 @@ void AliEMCALRecPoint::EvalLocalPosition(Float_t logWeight, TClonesArray * digit
   // Calculates the center of gravity in the local EMCAL-module coordinates 
   //  Info("Print", " logWeight %f : cluster energy %f ", logWeight, fAmp); // for testing
   
+  static Double_t dist;
+
   AliEMCALDigit * digit;
-  Int_t i=0, nstat=0;
+  Int_t i=0, nstat=0, idMax=-1;
   Double_t clXYZ[3]={0.,0.,0.}, clRmsXYZ[3]={0.,0.,0.}, xyzi[3], wtot=0., w=0.;
 
+  //printf(" dist : %f  e : %f \n", dist, fAmp);
   for(Int_t iDigit=0; iDigit<fMulDigit; iDigit++) {
     digit = dynamic_cast<AliEMCALDigit *>(digits->At(fDigitsList[iDigit])) ;
+    if(iDigit==0) {
+      idMax = digit->GetId(); // is it correct
+      dist  = TmaxInCm(Double_t(fAmp));
+    }
+    fGeomPtr->RelPosCellInSModule(digit->GetId(), idMax, dist, xyzi[0], xyzi[1], xyzi[2]);
+    //printf(" Id %i : dist %f Local x,y,z %f %f %f \n", digit->GetId(), dist, xyzi[0], xyzi[1], xyzi[2]);
 
-    fGeomPtr->RelPosCellInSModule(digit->GetId(), xyzi[0], xyzi[1], xyzi[2]);
-    // printf(" Id %i : Local x,y,z %f %f %f \n", digit->GetId(), xyzi[0], xyzi[1], xyzi[2]);
+    //fGeomPtr->RelPosCellInSModule(digit->GetId(), xyzi[0], xyzi[1], xyzi[2]);
+    //printf(" Id %i : dist %f Local x,y,z %f %f %f \n", digit->GetId(), 0.0, xyzi[0], xyzi[1], xyzi[2]);
+    //    if(fAmp>102.) assert(0);
 
     if(logWeight > 0.0)  w = TMath::Max( 0., logWeight + TMath::Log( fEnergyList[iDigit] / fAmp ));
     else  w = fEnergyList[iDigit]; // just energy
@@ -523,6 +536,163 @@ void AliEMCALRecPoint::EvalLocalPosition(Float_t logWeight, TClonesArray * digit
 //void AliEMCALRecPoint::EvalLocalPositionSimple()
 //{ // Weight is proportional of cell energy 
 //}
+//____________________________________________________________________________
+void AliEMCALRecPoint::EvalLocalPositionFit(Double_t deff, Double_t logWeight, 
+Double_t phiSlope, TClonesArray * digits)
+{
+  // Aug 14-16, 2007 - for fit 
+  // Aug 31 - should be static ??
+  static Double_t dist, ycorr;
+  static AliEMCALDigit *digit;
+
+  Int_t i=0, nstat=0, idMax=-1;
+  Double_t clXYZ[3]={0.,0.,0.}, clRmsXYZ[3]={0.,0.,0.}, xyzi[3], wtot=0., w=0.; 
+
+  for(Int_t iDigit=0; iDigit<digits->GetEntries(); iDigit++) {
+    digit = dynamic_cast<AliEMCALDigit *>(digits->At(fDigitsList[iDigit])) ;
+    if(iDigit==0) {
+      idMax = digit->GetId(); // is it correct
+      dist  = TmaxInCm(Double_t(fAmp));
+    }
+
+    dist = deff;
+    fGeomPtr->RelPosCellInSModule(digit->GetId(), idMax, dist, xyzi[0], xyzi[1], xyzi[2]);
+
+    if(logWeight > 0.0)  w = TMath::Max( 0., logWeight + TMath::Log( fEnergyList[iDigit] / fAmp ));
+    else                 w = fEnergyList[iDigit]; // just energy
+
+    if(w>0.0) {
+      wtot += w ;
+      nstat++;
+      for(i=0; i<3; i++ ) {
+        clXYZ[i]    += (w*xyzi[i]);
+        clRmsXYZ[i] += (w*xyzi[i]*xyzi[i]);
+      }
+    }
+  }
+  //  cout << " wtot " << wtot << endl;
+  if ( wtot > 0 ) { 
+    //    xRMS   = TMath::Sqrt(x2m - xMean*xMean);
+    for(i=0; i<3; i++ ) {
+      clXYZ[i] /= wtot;
+      if(nstat>1) {
+        clRmsXYZ[i] /= (wtot*wtot);  
+        clRmsXYZ[i]  =  clRmsXYZ[i] - clXYZ[i]*clXYZ[i];
+	if(clRmsXYZ[i] > 0.0) {
+          clRmsXYZ[i] =  TMath::Sqrt(clRmsXYZ[i]);
+        } else      clRmsXYZ[i] = 0;
+      } else        clRmsXYZ[i] = 0;
+    }    
+  } else {
+    for(i=0; i<3; i++ ) {
+      clXYZ[i] = clRmsXYZ[i] = -1.;
+    }
+  }
+  // clRmsXYZ[i] ??
+  if(phiSlope != 0.0 && logWeight > 0.0 && wtot) { 
+    // Correction in phi direction (y - coords here); Aug 16;
+    // May be put to global level or seperate method
+    ycorr = clXYZ[1] * (1. + phiSlope);
+    //printf(" y %f : ycorr %f : slope %f \n", clXYZ[1], ycorr, phiSlope); 
+    clXYZ[1] = ycorr;
+  }
+  fLocPos.SetX(clXYZ[0]);
+  fLocPos.SetY(clXYZ[1]);
+  fLocPos.SetZ(clXYZ[2]);
+    
+//  if (gDebug==2)
+//    printf("EvalLocalPosition: eta,phi,r = %f,%f,%f", fLocPos.X(), fLocPos.Y(), fLocPos.Z()) ; 
+  fLocPosM = 0 ; // covariance matrix
+}
+
+Bool_t AliEMCALRecPoint::EvalLocalPosition2(TClonesArray * digits, TArrayD &ed)
+{
+  // Evaluated local position of rec.point using digits 
+  // and parametrisation of w0 and deff
+  //printf(" <I> AliEMCALRecPoint::EvalLocalPosition2() \n"); 
+  return AliEMCALRecPoint::EvalLocalPositionFromDigits(digits, ed, fLocPos);
+}
+
+Bool_t AliEMCALRecPoint::EvalLocalPositionFromDigits(TClonesArray *digits, TArrayD &ed, TVector3 &locPos)
+{
+  // Used when digits should be recalibrated
+  static Double_t deff, w0, esum;
+  static Int_t iDigit;
+  static AliEMCALDigit *digit;
+
+  if(ed.GetSize() && (digits->GetEntries()!=ed.GetSize())) return kFALSE;
+
+  // Calculate sum energy of digits
+  esum = 0.0;
+  for(iDigit=0; iDigit<ed.GetSize(); iDigit++) esum += ed[iDigit];
+
+  GetDeffW0(esum, deff, w0);
+  
+  return EvalLocalPositionFromDigits(esum, deff, w0, digits, ed, locPos); 
+}
+
+Bool_t AliEMCALRecPoint::EvalLocalPositionFromDigits(const Double_t esum, const Double_t deff, const Double_t w0, 
+						     TClonesArray *digits, TArrayD &ed, TVector3 &locPos)
+{
+  static AliEMCALDigit *digit;
+
+  Int_t i=0, nstat=0, idMax=-1;
+  Double_t clXYZ[3]={0.,0.,0.}, xyzi[3], wtot=0., w=0.; 
+
+  AliEMCALGeometry* geo = AliEMCALGeometry::GetInstance(); // Get pointer to EMCAL geometry
+   
+  for(Int_t iDigit=0; iDigit<digits->GetEntries(); iDigit++) {
+    digit = dynamic_cast<AliEMCALDigit *>(digits->At(iDigit));
+
+    geo->RelPosCellInSModule(digit->GetId(), idMax, deff, xyzi[0], xyzi[1], xyzi[2]);
+
+    if(w0 > 0.0)  w = TMath::Max( 0., w0 + TMath::Log(ed[iDigit] / esum));
+    else          w = ed[iDigit]; // just energy
+
+    if(w>0.0) {
+      wtot += w ;
+      nstat++;
+      for(i=0; i<3; i++ ) {
+        clXYZ[i] += (w*xyzi[i]);
+      }
+    }
+  }
+  //  cout << " wtot " << wtot << endl;
+  if (wtot > 0) { 
+    for(i=0; i<3; i++ ) {
+      clXYZ[i] /= wtot;
+    }
+    locPos.SetX(clXYZ[0]);
+    locPos.SetY(clXYZ[1]);
+    locPos.SetZ(clXYZ[2]);
+    return kTRUE;
+  } else {
+    return kFALSE;
+  }
+
+}
+
+void AliEMCALRecPoint::GetDeffW0(const Double_t esum , Double_t &deff,  Double_t &w0)
+{
+  //
+  // Aug 31, 2001 
+  // Applied for simulation data with threshold 3 adc
+  // Calculate efective distance (deff) and weigh parameter (w0) 
+  // for coordinate calculation; 0.5 GeV < esum <100 GeV.
+  // Look to:  http://rhic.physics.wayne.edu/~pavlinov/ALICE/SHISHKEBAB/RES/CALIB/GEOMCORR/deffandW0VaEgamma_2.gif
+  //
+  static Double_t e=0.0;
+  const  Double_t dp0=9.25147, dp1=1.16700; // Hard coded now
+  const  Double_t wp0=4.83713, wp1=-2.77970e-01, wp2 = 4.41116;
+
+  // No extrapolation here
+  e = esum<0.5?0.5:esum;
+  e = e>100.?100.:e;
+
+  deff = dp0 + dp1*TMath::Log(e);
+  w0   = wp0 / (1. + TMath::Exp(wp1*(e+wp2)));
+  //printf("<I> AliEMCALRecPoint::GetDeffW0 esum %5.2f : deff %5.2f : w0 %5.2f \n", esum, deff, w0); 
+}
 
 //______________________________________________________________________________
 void AliEMCALRecPoint::EvalCoreEnergy(Float_t logWeight, TClonesArray * digits)
@@ -904,6 +1074,25 @@ void AliEMCALRecPoint::Paint(Option_t *)
   gPad->PaintPolyMarker(1,&x,&y,"") ;
 }
 
+Double_t AliEMCALRecPoint::TmaxInCm(const Double_t e , const Int_t key)
+{ 
+  // e energ in in GeV)
+  // key  =  0(gamma, default)
+  //     !=  0(electron)
+  static Double_t ca = 4.82;  // shower max parameter - first guess; ca=TMath::Log(1000./8.07)
+  static Double_t X0 = 1.23;  // radiation lenght (cm)
+  static Double_t tmax = 0.;   // position of electromagnetic shower max in cm
+
+  tmax = 0.0;
+  if(e>0.1) {
+    tmax = TMath::Log(e) + ca;
+    if      (key==0) tmax += 0.5; 
+    else             tmax -= 0.5;
+    tmax *= X0; // convert to cm
+  }
+  return tmax;
+}
+
 //______________________________________________________________________________
 Float_t AliEMCALRecPoint::EtaToTheta(Float_t arg) const
 {
diff --git a/EMCAL/AliEMCALRecPoint.h b/EMCAL/AliEMCALRecPoint.h
index f7d5f78f6d7..e2c53ce8b87 100644
--- a/EMCAL/AliEMCALRecPoint.h
+++ b/EMCAL/AliEMCALRecPoint.h
@@ -42,11 +42,19 @@ class AliEMCALRecPoint : public AliRecPoint {
   virtual Int_t   GetClusterType() const { return fClusterType; }
 
   virtual void    EvalAll(Float_t logWeight, TClonesArray * digits);
-  virtual void    EvalLocalPosition(Float_t logWeight, TClonesArray * digits) ;
-  //  void            EvalLocalPositionSimple(TClonesArray *digits); // ??
+  virtual void    EvalLocalPosition(Float_t logWeight, TClonesArray * digits);
   virtual void    EvalPrimaries(TClonesArray * digits) ;
   virtual void    EvalParents(TClonesArray * digits) ;
 
+  //  void            EvalLocalPositionSimple(TClonesArray *digits); // ??
+  void   EvalLocalPositionFit(Double_t deff, Double_t w0, Double_t phiSlope,TClonesArray * digits);
+  Bool_t EvalLocalPosition2(TClonesArray *digits, TArrayD &ed);
+
+  static  Bool_t  EvalLocalPositionFromDigits(const Double_t esum, const Double_t deff, const Double_t w0, 
+                                              TClonesArray *digits, TArrayD &ed, TVector3 &locPos);
+  static  Bool_t  EvalLocalPositionFromDigits(TClonesArray *digits, TArrayD &ed, TVector3 &locPos);
+  static  void    GetDeffW0(const Double_t esum, Double_t &deff,  Double_t &w0);
+
   using AliRecPoint::GetGlobalPosition;
   virtual void    GetGlobalPosition(TVector3 & gpos) const; // return global position (x, y, z) in ALICE
   virtual void    GetLocalPosition(TVector3 & lpos) const;  // return local position  (x, y, z) in EMCAL SM
@@ -67,6 +75,9 @@ class AliEMCALRecPoint : public AliRecPoint {
   Int_t       GetMultiplicityAtLevel(Float_t level) const ;  // computes multiplicity of digits with 
   Int_t *     GetAbsId() const {return fAbsIdList;}
   Int_t       GetAbsId(int i) const {if(i>=0 && i<fMulDigit)return fAbsIdList[i]; else return -1;}
+  Int_t       GetAbsIdMaxDigit() {return GetAbsId(fDigitIndMax);}
+  Int_t       GetIndMaxDigit() {return fDigitIndMax;}
+  void        SetIndMaxDigit(const Int_t ind) {fDigitIndMax = ind;}
   // energy above relative level
   virtual Int_t GetNumberOfLocalMax(AliEMCALDigit **  maxAt, Float_t * maxAtEnergy,
                                     Float_t locMaxCut,TClonesArray * digits ) const ; 
@@ -87,6 +98,7 @@ class AliEMCALRecPoint : public AliRecPoint {
     Fatal("operator =", "not implemented") ;
     return *this ; 
   }
+  static Double_t TmaxInCm(const Double_t e=0.0, const Int_t key=0);
 
 protected:
           void  EvalCoreEnergy(Float_t logWeight,TClonesArray * digits) ;             
@@ -104,7 +116,7 @@ private:
 				   // pseudocluster or v1
 	  Float_t fCoreEnergy ;       // energy in a shower core 
 	  Float_t fLambda[2] ;        // shower ellipse axes
-	  Float_t fDispersion ;       // shower dispersio
+	  Float_t fDispersion ;       // shower dispersion
 	  Float_t *fEnergyList ;      //[fMulDigit] energy of digits
 	  Float_t *fTimeList ;        //[fMulDigit] time of digits
           Int_t   *fAbsIdList;        //[fMulDigit] absId  of digits
@@ -116,8 +128,9 @@ private:
           Int_t * fParentsList;       // [fMulParent] list of the parents of the digits
           Float_t * fDEParentsList;   // [fMulParent] list of the parents of the digits
           Int_t   fSuperModuleNumber; // number identifying supermodule containing recpoint
-
-  ClassDef(AliEMCALRecPoint,8) // RecPoint for EMCAL (Base Class)
+          // Aug 16, 2007
+          Int_t   fDigitIndMax;       // Index of digit with max energy in array fAbsIdList
+  ClassDef(AliEMCALRecPoint,9) // RecPoint for EMCAL (Base Class)
  
 };
 
diff --git a/EMCAL/AliEMCALShishKebabTrd1Module.cxx b/EMCAL/AliEMCALShishKebabTrd1Module.cxx
index ab581048b3e..678a0ca116c 100644
--- a/EMCAL/AliEMCALShishKebabTrd1Module.cxx
+++ b/EMCAL/AliEMCALShishKebabTrd1Module.cxx
@@ -21,7 +21,7 @@
 // Sep 20004 - Nov 2006
 // See web page with description of Shish-Kebab geometries:
 // http://pdsfweb01.nersc.gov/~pavlinov/ALICE/SHISHKEBAB/RES/shishkebabALICE.html
-// Nov 9,2006 - added cas of 3X3
+// Nov 9,2006 - added case of 3X3
 //_________________________________________________________________________
 
 #include "AliLog.h"
@@ -54,6 +54,8 @@ AliEMCALShishKebabTrd1Module::AliEMCALShishKebabTrd1Module(Double_t theta, AliEM
     fOB(),
     fOB1(),
     fOB2(),
+    fThetaOB1(0.),
+    fThetaOB2(0.),
     fOK3X3()
 { 
   // theta in radians ; first object shold be with theta=pi/2.
@@ -81,6 +83,8 @@ AliEMCALShishKebabTrd1Module::AliEMCALShishKebabTrd1Module(AliEMCALShishKebabTrd
     fOB(),
     fOB1(),
     fOB2(),
+    fThetaOB1(0.),
+    fThetaOB2(0.),
     fOK3X3()
 { 
   //  printf("** Left Neighbor : %s **\n", leftNeighbor.GetName());
@@ -141,7 +145,7 @@ void AliEMCALShishKebabTrd1Module::Init(Double_t A, Double_t B)
 void AliEMCALShishKebabTrd1Module::DefineAllStaff()
 {
   DefineName(fTheta);
-  // Centers of module - 2X2 case
+  // Centers of cells - 2X2 case
   Double_t kk1 = (fga+fga2)/(2.*4.); // kk1=kk2 
 
   Double_t xk1 = fOK.X() - kk1*TMath::Sin(fTheta);
@@ -152,7 +156,7 @@ void AliEMCALShishKebabTrd1Module::DefineAllStaff()
   Double_t yk2 = fOK.Y() - kk1*TMath::Cos(fTheta) - fgr;
   fOK2.Set(xk2,yk2);
 
-  // Centers of module - 3X3 case; Nov 9,2006
+  // Centers of cells - 3X3 case; Nov 9,2006
   fOK3X3[1].Set(fOK.X(), fOK.Y()-fgr); // coincide with module center
 
   kk1 = ((fga+fga2)/4. + fga/6.)/2.; 
@@ -165,11 +169,13 @@ void AliEMCALShishKebabTrd1Module::DefineAllStaff()
   yk2 = fOK.Y() - kk1*TMath::Cos(fTheta) - fgr;
   fOK3X3[2].Set(xk2,yk2);
 
-  // May 15, 2006; position of cell face of cells 
+  // May 15, 2006; position of module(cells) center face 
   fOB.Set(fOK.X()-fgb/2.*TMath::Cos(fTheta),  fOK.Y()-fgb/2.*TMath::Sin(fTheta)-fgr);
   fOB1.Set(fOB.X()-fga/4.*TMath::Sin(fTheta), fOB.Y()+fga/4.*TMath::Cos(fTheta));
   fOB2.Set(fOB.X()+fga/4.*TMath::Sin(fTheta), fOB.Y()-fga/4.*TMath::Cos(fTheta));
-
+  // Jul 30, 2007 - for taking into account a position of shower maximum
+  fThetaOB1 = fTheta - fgangle/4.; // ??
+  fThetaOB2 = fTheta + fgangle/4.;
 }
 
 //_____________________________________________________________________________
@@ -264,6 +270,26 @@ Double_t  AliEMCALShishKebabTrd1Module::GetEtaOfCenterOfModule() const
   return -TMath::Log(TMath::Tan(fOK.Phi()/2.));
 }
 
+void AliEMCALShishKebabTrd1Module::GetPositionAtCenterCellLine(Int_t ieta, Double_t dist, TVector2 &v)
+{
+  // Jul 30, 2007
+  static Double_t theta=0., x=0., y=0.;
+  if(ieta==0) {
+    v     = fOB2;
+    theta = fThetaOB2;
+  } else if(ieta==1) {
+    v     = fOB1;
+    theta = fThetaOB1;
+  } else {
+    assert(0);
+  }
+  x = v.X() + TMath::Cos(theta) * dist;
+  y = v.Y() + TMath::Sin(theta) * dist;
+  v.Set(x,y);
+  //printf(" GetPositionAtCenterCellLine() : dist %f : ieta %i : x %f | y %f \n", dist, ieta, x, y);
+}  
+
+
 //_____________________________________________________________________________
 Double_t  AliEMCALShishKebabTrd1Module::GetMaxEtaOfModule(int pri) const 
 { 
diff --git a/EMCAL/AliEMCALShishKebabTrd1Module.h b/EMCAL/AliEMCALShishKebabTrd1Module.h
index 730a4cc440e..859ae44ba77 100644
--- a/EMCAL/AliEMCALShishKebabTrd1Module.h
+++ b/EMCAL/AliEMCALShishKebabTrd1Module.h
@@ -73,6 +73,8 @@ class AliEMCALShishKebabTrd1Module : public TNamed {
     if(ieta<=0) return fOB2;
     else        return fOB1;
   }  
+  // Jul 30, 2007
+  void GetPositionAtCenterCellLine(Int_t ieta, Double_t dist, TVector2 &v);
   // 
   Double_t GetTanBetta() const {return fgtanBetta;}
   Double_t Getb()        const {return fgb;}
@@ -103,14 +105,17 @@ class AliEMCALShishKebabTrd1Module : public TNamed {
   // Nov 04,2004; Feb 19,2006 
   TVector2 fOK1; // ieta=1
   TVector2 fOK2; // ieta=0
-  // May 13, 2006; local position of module (cells) face  
-  TVector2 fOB;  // module
-  TVector2 fOB1; // ieta=1
-  TVector2 fOB2; // ieta=0
+  // May 13, 2006; local position of module (cells) center face  
+  TVector2 fOB;   // module
+  TVector2 fOB1;  // ieta=1
+  TVector2 fOB2;  // ieta=0
+  // Jul 30, 2007
+  Double_t fThetaOB1; // theta of cell center line (go through OB1)
+  Double_t fThetaOB2;  // theta of cell center line (go through OB2)
   // 3X3 case - Nov 9,2006
   TVector2 fOK3X3[3];
   // public:
-  ClassDef(AliEMCALShishKebabTrd1Module,1) // TRD1 Shish-Kebab module 
+  ClassDef(AliEMCALShishKebabTrd1Module,2) // TRD1 Shish-Kebab module 
 };
 
 #endif
diff --git a/EMCAL/AliEMCALv0.cxx b/EMCAL/AliEMCALv0.cxx
index 599ab399d17..ef236318a7b 100644
--- a/EMCAL/AliEMCALv0.cxx
+++ b/EMCAL/AliEMCALv0.cxx
@@ -743,7 +743,8 @@ void AliEMCALv0::CreateEmod(const char* mother, const char* child)
       fParEMOD[2] = g->GetPhiModuleSize()/2.;;  // dy
       fParEMOD[3] = g->GetLongModuleSize()/2.;  // dz
       gMC->Gsvolu(child, "TRD1", fIdTmedArr[kIdSTEEL], fParEMOD, 4);
-      if(gn.Contains("WSUC") || gn.Contains("MAY05")){
+      //      if(gn.Contains("WSUC") || gn.Contains("MAY05")){
+      if(0){ // Jul 12 - should be checked
         parSCPA[0] = g->GetEtaModuleSize()/2. + tanTrd1*g->GetFrontSteelStrip();   // dx1
         parSCPA[1] = parSCPA[0]               + tanTrd1*g->GetPassiveScintThick(); // dx2
         parSCPA[2] = g->GetPhiModuleSize()/2.;     // dy
@@ -1299,6 +1300,17 @@ AliEMCALShishKebabTrd1Module* AliEMCALv0::GetShishKebabModule(Int_t neta)
 
 //_____________________________________________________________________________
 void AliEMCALv0::AddAlignableVolumes() const
+{
+  TString ntmp(GetTitle()); // name of EMCAL geometry
+
+  if(ntmp.Contains("WSUC")) {
+    AddAlignableVolumesInWSUC();  // WSUC case
+  } else {
+    AddAlignableVolumesInALICE(); // ALICE case
+  }
+}
+
+void AliEMCALv0::AddAlignableVolumesInALICE() const
 {
   //
   // Create entries for alignable volumes associating the symbolic volume
@@ -1337,3 +1349,25 @@ void AliEMCALv0::AddAlignableVolumes() const
 
 }
 
+void AliEMCALv0::AddAlignableVolumesInWSUC() const
+{
+  //
+  // Create entries for alignable volumes associating the symbolic volume
+  // name with the corresponding volume path. Needs to be synchronized with
+  // eventual changes in the geometry.
+  //
+
+  TString vpstr1 = "WSUC_1/XEN1_1/SMOD_";
+  TString snstr1 = "EMCAL/CosmicTestSupermodule";
+  TString volpath, symname;
+
+  // #SM is just one 
+  for (Int_t smodnum=0; smodnum < 1; smodnum++) {
+    symname = snstr1;
+    symname += (smodnum+1);
+    volpath = vpstr1;
+    volpath += (smodnum+1);
+    if(!gGeoManager->SetAlignableEntry(symname.Data(),volpath.Data()))
+      AliFatal("Unable to set alignable entry!!");
+  }
+}
diff --git a/EMCAL/AliEMCALv0.h b/EMCAL/AliEMCALv0.h
index 0fd2239954a..705fd867e76 100644
--- a/EMCAL/AliEMCALv0.h
+++ b/EMCAL/AliEMCALv0.h
@@ -36,6 +36,8 @@ class AliEMCALv0 : public AliEMCAL {
   using AliEMCAL::AddHit;
 
   virtual void  AddAlignableVolumes() const;
+  virtual void  AddAlignableVolumesInALICE() const;
+  virtual void  AddAlignableVolumesInWSUC() const;
 
   virtual void BuildGeometry();// creates the geometry for the ROOT display
   TNode *BuildGeometryOfWSUC();  // WSUC - test environment
diff --git a/EMCAL/EMCALLinkDefbase.h b/EMCAL/EMCALLinkDefbase.h
index 5c93d52b5dd..f271546c2bf 100644
--- a/EMCAL/EMCALLinkDefbase.h
+++ b/EMCAL/EMCALLinkDefbase.h
@@ -25,14 +25,17 @@
 #pragma link C++ class AliEMCALCalibHistoProducer+;
 #pragma link C++ class AliEMCALPreprocessor+;
 #pragma link C++ class AliEMCALRawUtils+;
+// Calibration staff
+#pragma link C++ class AliEMCALFolder;
+#pragma link C++ class AliEMCALSuperModule;
+#pragma link C++ class AliEMCALCell;
+// Tables without TTable
+#pragma link C++ class AliEMCALPi0SelectionParam+;
+#pragma link C++ class pi0SelectionParam+;
+#pragma link C++ class AliEMCALCalibCoefs+;
+#pragma link C++ class calibCoef+;
+#pragma link C++ class AliEMCALCellInfo+;
+#pragma link C++ class cellInfo+;
+//
 #pragma link C++ class AliCaloCalibPedestal+;
-// Calibration stuff
-// #pragma link C++ class AliEMCALFolder;
-// #pragma link C++ class AliEMCALSuperModule;
-// #pragma link C++ class AliEMCALCell;
-// Tables
-// #pragma link C++ class AliEMCALPi0SelectionParam-;
-// #pragma link C++ class pi0SelectionParam+;
-// #pragma link C++ class AliEMCALCalibCoefs-;
-// #pragma link C++ class calibCoef+;
 #endif
diff --git a/EMCAL/EMCALLinkDefrec.h b/EMCAL/EMCALLinkDefrec.h
index 076e699e8aa..1b9fe5fb8ed 100644
--- a/EMCAL/EMCALLinkDefrec.h
+++ b/EMCAL/EMCALLinkDefrec.h
@@ -9,7 +9,8 @@
 #pragma link C++ class AliEMCALTrack+;
 #pragma link C++ class AliEMCALTracker+;
 #pragma link C++ class AliEMCALPID+;
-// #pragma link C++ class AliEMCALRecPointsQaESDSelector+;
+#pragma link C++ class AliEMCALRecPointsQaESDSelector+;
+//
 #pragma link C++ class AliEMCALRecParam+;
 
 #endif
diff --git a/EMCAL/libEMCALbase.pkg b/EMCAL/libEMCALbase.pkg
index a416883f41e..9c3a672911b 100644
--- a/EMCAL/libEMCALbase.pkg
+++ b/EMCAL/libEMCALbase.pkg
@@ -20,12 +20,13 @@ AliEMCALCalibData.cxx \
 AliEMCALCalibHistoProducer.cxx \
 AliEMCALPreprocessor.cxx \
 AliEMCALRawUtils.cxx \
-AliCaloCalibPedestal.cxx \
-# AliEMCALSuperModule.cxx \
-# AliEMCALCell.cxx \
-# AliEMCALFolder.cxx \
-# AliEMCALPi0SelectionParam.cxx \
-# AliEMCALCalibCoefs.cxx
+AliEMCALFolder.cxx \
+AliEMCALSuperModule.cxx \
+AliEMCALCell.cxx \
+AliEMCALPi0SelectionParam.cxx \
+AliEMCALCalibCoefs.cxx \
+AliEMCALCellInfo.cxx\
+AliCaloCalibPedestal.cxx
 
 HDRS= $(SRCS:.cxx=.h) 
 
diff --git a/EMCAL/libEMCALrec.pkg b/EMCAL/libEMCALrec.pkg
index 7998af4e48c..42dd0bd2772 100644
--- a/EMCAL/libEMCALrec.pkg
+++ b/EMCAL/libEMCALrec.pkg
@@ -6,7 +6,7 @@ AliEMCALTrack.cxx \
 AliEMCALTracker.cxx \
 AliEMCALPID.cxx \
 AliEMCALRecParam.cxx \
-# AliEMCALRecPointsQaESDSelector.cxx \
+AliEMCALRecPointsQaESDSelector.cxx \
 
 HDRS= $(SRCS:.cxx=.h)