HLT/EMCAL/AliHLTEMCALGeometry.cxx

   1 /**************************************************************************
   2  * This file is property of and copyright by the Experimental Nuclear     *
   3  * Physics Group, Dep. of Physics                                         *
   4  * University of Oslo, Norway, 2007                                       *
   5  *                                                                        *
   6  * Author: federico ronchetti         for the ALICE HLT Project.*
   7  *                                                                        *
   8  * Permission to use, copy, modify and distribute this software and its   *
   9  * documentation strictly for non-commercial purposes is hereby granted   *
  10  * without fee, provided that the above copyright notice appears in all   *
  11  * copies and that both the copyright notice and this permission notice   *
  12  * appear in the supporting documentation. The authors make no claims     *
  13  * about the suitability of this software for any purpose. It is          *
  14  * provided "as is" without express or implied warranty.                  *
  15  **************************************************************************/
  16
  17 #include "AliHLTEMCALGeometry.h"
  18 #include "AliHLTEMCALConstants.h"
  19 #include "AliCDBEntry.h"
  20 #include "AliCDBManager.h"
  21 #include "AliCDBPath.h"
  22
  23
  24 ClassImp(AliHLTEMCALGeometry);
  25 TGeoManager *gGeoManager = 0;
  26
  27 AliHLTEMCALGeometry::AliHLTEMCALGeometry() :
  28         AliHLTCaloGeometry ("EMCAL"),
  29         fGeo(0)
  30 {
  31
  32   //fGeo = new AliEMCALGeoUtils("EMCAL_COMPLETE","EMCAL");
  33         //fGeo = new AliEMCALGeometry("EMCAL_COMPLETE","EMCAL");
  34         //fGeo =  AliEMCALGeometry::GetInstance(AliEMCALGeometry::GetDefaultGeometryName());
  35         //TGeoManager::Import("/home/fedro/work/AliRoot/test/QA/geometry.root");
  36
  37   GetGeometryFromCDB();
  38 }
  39
  40 Int_t AliHLTEMCALGeometry::InitialiseGeometry()
  41 {
  42
  43    return GetGeometryFromCDB();
  44 }
  45
  46
  47 AliHLTEMCALGeometry::~AliHLTEMCALGeometry()
  48 {
  49
  50 }
  51
  52 void
  53 AliHLTEMCALGeometry::GetGlobalCoordinates(AliHLTCaloRecPointDataStruct &recPoint, AliHLTCaloGlobalCoordinate &globalCoord)
  54 {
  55         // Working variables
  56   Int_t istrip = 0;
  57   Float_t z0 = 0;
  58   Float_t zb = 0;
  59   Float_t z_is = 0;
  60   Float_t d = 0;
  61   Float_t x,y,z; // return variables in terry's RF
  62   Float_t dz = fCaloConstants->GetMINCELLSTEPETA(); // base cell width in eta
  63   Float_t dx = fCaloConstants->GetCELLSTEPPHI(); // base cell width in phi
  64
  65   // parameters for shower depth calculation
  66   Float_t X0 = fCaloConstants->GetRADLENGTH();
  67   Float_t Ecr = fCaloConstants->GetCRITICENERGY();
  68   Float_t Cj = fCaloConstants->GetCJ();
  69
  70   // tapering angles
  71   Float_t teta0 = fCaloConstants->GetCELLANGLE(); //tapering angle (deg)
  72   Float_t teta1; //working angle
  73   Float_t L = fCaloConstants->GetCELLHEIGHT();
  74
  75   // converting to MEV
  76   Float_t E = recPoint.fAmp * 1000;
  77
  78   Double_t glob[3] ={0,0,0};
  79   Double_t loc[3] = {0,0,0};
  80
  81   if (recPoint.fZ >= 47.51 || recPoint.fZ< -0.51) {
  82     //Logging(kHLTLogError, "HLT", "EMCAL", "AliHLTEMCALGeometry::GetGlobalCoordinates: invalid Z: %f from recpoint ", recPoint.fZ);
  83     return;
  84   }
  85
  86   if (recPoint.fX >= 23.51 || recPoint.fX< -0.51) {
  87     //Logging(kHLTLogError, "HLT", "EMCAL", "AliHLTEMCALGeometry::GetGlobalCoordinates: invalid X: % from recpoint ", recPoint.fX);
  88     return;
  89   }
  90
  91
  92   switch ( recPoint.fParticle )
  93     {
  94     case 0:
  95       Cj = - fCaloConstants->GetCJ(); // photon
  96       d = X0 * TMath::Log( E / Ecr + Cj);
  97       break;
  98
  99     case 1:
 100       Cj = + fCaloConstants->GetCJ(); // electron
 101       d = X0 * TMath::Log( E / Ecr + Cj);
 102       break;
 103
 104     case 2:
 105       // hadron
 106       d = 0.5 * L;
 107       break;
 108
 109     default:
 110       Cj = - fCaloConstants->GetCJ(); // defaulting to photon
 111       d = X0 * TMath::Log( E / Ecr + Cj);
 112     }
 113
 114   istrip = int ((recPoint.fZ + 0.5 ) / 2);
 115
 116   // module angle
 117   teta1 = TMath::DegToRad() * istrip * teta0;
 118
 119   // calculation of module corner along z
 120   // as a function of strip
 121
 122   for (Int_t is=0; is<= istrip; is++) {
 123
 124     teta1 = TMath::DegToRad() * is * teta0;
 125
 126     z_is = z_is + 2*dz*(TMath::Sin(teta1)*TMath::Tan(teta1) + TMath::Cos(teta1));
 127
 128   }
 129
 130   z0 = dz * (recPoint.fZ - 2*istrip + 0.5);
 131   zb = (2*dz-z0-d*TMath::Tan(teta1))*TMath::Cos(teta1);
 132
 133   z = z_is - zb*TMath::Cos(teta1);
 134
 135   //       cout << "----> istrip: " << istrip << endl;
 136   //       cout << "----> z0: "<< z0 << endl;
 137   //       cout << "----> zb: "<< zb << endl;
 138   //       cout << "----> corner z: "<< z_is << endl;
 139   //       cout << "----> teta1: "<< TMath::RadToDeg()*teta1 << endl;
 140
 141   y = d/TMath::Cos(teta1) + zb*TMath::Sin(teta1);
 142
 143   x = (recPoint.fX + 0.5)*dx;
 144
 145   // cout << "x: " << x << " y: "<< y << " z " << z << endl;
 146
 147   // check coordinate origin
 148   loc[0] = x;
 149   loc[1] = y;
 150   loc[2] = z;
 151
 152   if(!fGeo)
 153     {
 154       Logging(kHLTLogError, "HLT", "EMCAL", "AliHLTEMCALGeometry::GetGlobalCoordinates: no geometry initialised");
 155       return;
 156     }
 157
 158   //cout << "recpoint Particle " << recPoint.fParticle << endl;
 159   // cout << "recpoint Module " << recPoint.fModule << endl;
 160
 161          fGeo->GetGlobal(loc, glob, recPoint.fModule);
 162
 163
 164          // cout << "recpoint Amp " << recPoint.fAmp << endl;
 165
 166          ConvertRecPointCoordinates(loc[0], loc[1], loc[2]);
 167
 168          //cout << "after conversion ";
 169          //cout << "x: " << loc[0] << " y: "<< loc[1] << " z " << loc[2] << endl;
 170
 171          fGeo->GetGlobal(&loc[0], &glob[0], recPoint.fModule);
 172
 173          //cout << "global cooordinates after call:";
 174          //cout << "x: " << glob[0] << " y: "<< glob[1] << " z " << glob[2] << endl;
 175
 176          globalCoord.fX = glob[0];
 177          globalCoord.fY = glob[1];
 178          globalCoord.fZ = glob[2];
 179 }
 180
 181 void
 182 AliHLTEMCALGeometry::GetCellAbsId(UInt_t module, UInt_t x, UInt_t z, Int_t& AbsId)
 183 {
 184
 185   if(!fGeo)
 186     {
 187       Logging(kHLTLogError, "HLT", "EMCAL", "AliHLTEMCALGeometry::GetCellAbsId: no geometry initialised");
 188       return;
 189
 190     }
 191         AbsId = fGeo->GetAbsCellIdFromCellIndexes(module, (Int_t) x, (Int_t) z);
 192
 193
 194
 195 }
 196
 197 void AliHLTEMCALGeometry::ConvertRecPointCoordinates(Double_t &x, Double_t &y, Double_t &z) const
 198 {
 199
 200         //FIXME
 201         // this should be read from GEANT box dimensions
 202         Double_t DX = 13.869008;
 203         Double_t DY = 72.559998;
 204         Double_t DZ = 175.00000;
 205
 206         // from our local frame to GEANT local frame
 207         x = y - DX;  //FIXME
 208         y = -x + DY; //FIXME
 209         z = z - DZ;  //FIXME
 210
 211 }
 212
 213
 214 int
 215 AliHLTEMCALGeometry::GetGeometryFromCDB()
 216 {
 217
 218   // AliCDBManager::Instance()->SetDefaultStorage("local://$ALICE_ROOT/OCDB");
 219
 220   AliCDBPath path("GRP","Geometry","Data");
 221   if(path.GetPath())
 222     {
 223       //      HLTInfo("configure from entry %s", path.GetPath());
 224       AliCDBEntry *pEntry = AliCDBManager::Instance()->Get(path/*,GetRunNo()*/);
 225       if (pEntry)
 226         {
 227           if(!fGeo)
 228             {
 229               delete fGeo;
 230               fGeo = 0;
 231             }
 232
 233           gGeoManager = (TGeoManager*) pEntry->GetObject();
 234
 235           if(gGeoManager)
 236             {
 237               fGeo = new AliEMCALGeoUtils("EMCAL_COMPLETE","EMCAL");
 238             }
 239
 240         }
 241       else
 242         {
 243           //HLTError("can not fetch object \"%s\" from OCDB", path);
 244           Logging(kHLTLogError, "HLT", "EMCAL", "can not fetch object from OCDB");
 245
 246         }
 247     }
 248   return 0;
 249 }