Tracking 2 Local system matrix added to GeoManager
authormarian <marian@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 22 May 2007 09:59:26 +0000 (09:59 +0000)
committermarian <marian@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 22 May 2007 09:59:26 +0000 (09:59 +0000)
Possible to use AliCluster::GetGlobalXYZ
The AliTPCDispalyClusterMI use GetGlobalXYZ  (Marian)

TPC/AliTPCDisplayClustersMI.C
TPC/AliTPCParam.cxx
TPC/AliTPCParam.h
TPC/AliTPCclusterMI.cxx
TPC/AliTPCclusterMI.h
TPC/AliTPCv2.cxx

index a3d10ec..432c05e 100644 (file)
@@ -1,7 +1,7 @@
 #if !defined(__CINT__) || defined(__MAKECINT__)
 #include "alles.h"
 #include "AliTPCtracker.h"
-#include "TView.h"
+#include "TView3D.h"
 #include "TPolyMarker3D.h"
 #include "AliSimDigits.h"
 #include "AliTPCParam.h"
@@ -13,6 +13,7 @@
   Author:   marian.ivanov@cern.ch
 
   How to use ?
+  TGeoManager::Import("geometry.root");
   .L AliTPCDisplayClustersMI.C+
   AliTPCDisplayClusters disp;
   disp.Init(0,0);   //specify event number, and threshold for the noise
@@ -123,7 +124,7 @@ void AliTPCDisplayClusters::DisplayClusters(Int_t first, Int_t last)
   Int_t MASK = 10000000;
 
   TCanvas *c1=new TCanvas("cdisplay", "Cluster display",0,0,700,730);
-  TView *v=new TView(1);
+  TView3D *v=new TView3D();
   v->SetRange(-330,-360,-330,360,360,1710);
   c1->Clear();
   c1->SetFillColor(1);
@@ -152,14 +153,17 @@ void AliTPCDisplayClusters::DisplayClusters(Int_t first, Int_t last)
        Float_t cs, sn, tmp;
        fParam->AdjustCosSin(sec,cs,sn);
        tmp = x*cs-y*sn; y= x*sn+y*cs; x=tmp; 
+       Float_t clxyz[3];
+       cl->GetGlobalXYZ(clxyz);
+
        Int_t trackId = cl->GetLabel(0);
        if ( (last>0) &&trackId>last) continue;
        if (trackId<first) continue;
        if (trackId < MASK-1) {
-         pmSignal->SetPoint(imarSignal,x,y,z);
+         pmSignal->SetPoint(imarSignal,clxyz[0],clxyz[1],clxyz[2]);
          imarSignal++;
        } else {
-         pm->SetPoint(imarBgr,x,y,z);
+         pm->SetPoint(imarBgr,clxyz[0],clxyz[1],clxyz[2]);
          imarBgr++;
        }          
       }
@@ -178,21 +182,21 @@ void AliTPCDisplayClusters::DisplayClusters(Int_t first, Int_t last)
   }
   
   
-  TNode * main = (TNode*)((fGeom->GetListOfNodes())->First());
-  TIter next(main->GetListOfNodes());
-  TNode  *module=0;
-  while((module = (TNode*)next())) {
-    char ch[100];
-    sprintf(ch,"%s\n",module->GetTitle());
-    //printf("%s\n",module->GetTitle());
-    if (ch[0]=='T'&&ch[1]=='P' && ch[2]=='C')  //if TPC draw
-      module->SetVisibility(3);
-    else
-      module->SetVisibility(-1);
-  }
+ //  TNode * main = (TNode*)((fGeom->GetListOfNodes())->First());
+//   TIter next(main->GetListOfNodes());
+//   TNode  *module=0;
+//   while((module = (TNode*)next())) {
+//     char ch[100];
+//     sprintf(ch,"%s\n",module->GetTitle());
+//     //printf("%s\n",module->GetTitle());
+//     if (ch[0]=='T'&&ch[1]=='P' && ch[2]=='C')  //if TPC draw
+//       module->SetVisibility(3);
+//     else
+//       module->SetVisibility(-1);
+//   }
   
   
-  fGeom->Draw("same");
+//   fGeom->Draw("same");
   c1->Modified(); c1->Update(); 
   
   
index b4f0982..e3e242c 100644 (file)
@@ -683,6 +683,19 @@ Bool_t AliTPCParam::ReadGeoMatrices(){
   return kTRUE;
 }
 
+TGeoHMatrix *  AliTPCParam::Tracking2LocalMatrix(const TGeoHMatrix * geoMatrix, Int_t sector) const{
+  //
+  // make local to tracking matrix
+  //
+  Double_t sectorAngle = 20.*(sector%18)+10;
+  TGeoHMatrix *newMatrix = new TGeoHMatrix();
+  newMatrix->RotateZ(sectorAngle);
+  newMatrix->MultiplyLeft(&(geoMatrix->Inverse()));
+  return newMatrix;
+}
+
+
+
 
 Bool_t AliTPCParam::GetStatus() const
 {
index 9bd6221..8243161 100644 (file)
@@ -24,7 +24,7 @@ class AliTPCParam : public AliDetectorParam {
 public:
   AliTPCParam(); 
   virtual ~AliTPCParam();
-  
+  TGeoHMatrix *  Tracking2LocalMatrix(const TGeoHMatrix * geoMatrix, Int_t sector) const;  
   virtual Bool_t  Transform(Float_t *xyz, Int_t *index, Int_t* oindex);
   //transformation from input coodination system to output coordination system  
   Int_t  Transform0to1(Float_t *xyz, Int_t *index) const;
index 49d8ba6..eba4a55 100644 (file)
@@ -28,6 +28,7 @@
 
 #include "AliTPCclusterMI.h"
 #include "AliTPCclusterInfo.h"
+#include "AliGeomManager.h"
 #include "AliLog.h"
 
 ClassImp(AliTPCclusterMI)
@@ -135,3 +136,15 @@ Int_t AliTPCclusterMI::Compare(const TObject* obj) const
   AliTPCclusterMI * o2 = (AliTPCclusterMI*)obj;
   return (o2->GetY()>GetY())? -1:1; 
 }
+
+
+void AliTPCclusterMI::SetDetector(Int_t detector){
+  //
+  // set volume ID 
+  //  
+  fDetector = (UChar_t)(detector%72);
+  AliGeomManager::ELayerID id = (fDetector<36) ? 
+    AliGeomManager::kTPC1 :AliGeomManager::kTPC2 ;
+  Int_t modId = (fDetector<36)?fDetector: fDetector-36;
+  SetVolumeId(AliGeomManager::LayerToVolUID(id,modId));  
+}
index cf776e9..aa975fc 100644 (file)
@@ -26,7 +26,7 @@ public:
   inline  void Use(Int_t inc=10);
   virtual Int_t GetDetector() const {return fDetector;}
   virtual Int_t GetRow() const {return fRow;}
-  virtual void SetDetector(Int_t detector){fDetector = (UChar_t)(detector%256);}
+  virtual void SetDetector(Int_t detector);
   virtual void SetRow(Int_t row){fRow = (UChar_t)(row%256);}  
   virtual void SetTimeBin(Float_t timeBin){ fTimeBin= timeBin;}
   virtual void SetPad(Float_t pad){ fPad = pad;}
index 6832f29..d5d1428 100644 (file)
@@ -45,6 +45,8 @@
 #include "TGeoTrd1.h"
 #include "TGeoCompositeShape.h"
 #include "TGeoPara.h"
+#include "TGeoPhysicalNode.h"
+
 ClassImp(AliTPCv2)
  
 //_____________________________________________________________________________
@@ -912,17 +914,25 @@ void AliTPCv2::SetInnerChambersAlignable() const
   TString snappend="/InnerChamber";
   TString volpath, symname;
   
-    for(Int_t cnt=1; cnt<=18; cnt++){
-      volpath = vpstr1;
-      volpath += cnt;
-      volpath += vpappend;
-      symname = snstr1;
-      symname += cnt;
-      symname += snappend;
+  for(Int_t cnt=1; cnt<=18; cnt++){
+    volpath = vpstr1;
+    volpath += cnt;
+    volpath += vpappend;
+    symname = snstr1;
+    symname += cnt;
+    symname += snappend;
     if(!gGeoManager->SetAlignableEntry(symname.Data(),volpath.Data()))
       AliFatal(Form("Alignable entry %s not created. Volume path %s not valid", symname.Data(),volpath.Data()));
-      modnum++;
-    }
+    //
+    TGeoPNEntry *alignableEntry = gGeoManager->GetAlignableEntry(symname.Data());
+    const char *path = alignableEntry->GetTitle();
+    if (!gGeoManager->cd(path))
+      AliFatal(Form("Volume path %s not valid!",path));
+    TGeoHMatrix* globMatrix = gGeoManager->GetCurrentMatrix();
+    TGeoHMatrix* matTtoL = fTPCParam->Tracking2LocalMatrix(globMatrix,cnt-1);
+    alignableEntry->SetMatrix(matTtoL);
+    modnum++;
+  }
 
   for(Int_t cnt=1; cnt<=18; cnt++){
     volpath = vpstr2;
@@ -933,6 +943,13 @@ void AliTPCv2::SetInnerChambersAlignable() const
     symname += snappend;
     if(!gGeoManager->SetAlignableEntry(symname.Data(),volpath.Data()))
       AliFatal(Form("Alignable entry %s not created. Volume path %s not valid", symname.Data(),volpath.Data()));
+    TGeoPNEntry *alignableEntry = gGeoManager->GetAlignableEntry(symname.Data());
+    const char *path = alignableEntry->GetTitle();
+    if (!gGeoManager->cd(path))
+      AliFatal(Form("Volume path %s not valid!",path));
+    TGeoHMatrix* globMatrix = gGeoManager->GetCurrentMatrix();
+    TGeoHMatrix* matTtoL = fTPCParam->Tracking2LocalMatrix(globMatrix,18+cnt-1);
+    alignableEntry->SetMatrix(matTtoL);
     modnum++;
   }
 }
@@ -959,6 +976,13 @@ void AliTPCv2::SetOuterChambersAlignable() const
     symname += snappend;
     if(!gGeoManager->SetAlignableEntry(symname.Data(),volpath.Data()))
       AliFatal(Form("Alignable entry %s not created. Volume path %s not valid", symname.Data(),volpath.Data()));
+    TGeoPNEntry *alignableEntry = gGeoManager->GetAlignableEntry(symname.Data());
+    const char *path = alignableEntry->GetTitle();
+    if (!gGeoManager->cd(path))
+      AliFatal(Form("Volume path %s not valid!",path));
+    TGeoHMatrix* globMatrix = gGeoManager->GetCurrentMatrix();
+    TGeoHMatrix* matTtoL = fTPCParam->Tracking2LocalMatrix(globMatrix,36+cnt-1);
+    alignableEntry->SetMatrix(matTtoL);
     modnum++;
   }
 
@@ -971,6 +995,13 @@ void AliTPCv2::SetOuterChambersAlignable() const
     symname += snappend;
     if(!gGeoManager->SetAlignableEntry(symname.Data(),volpath.Data()))
       AliFatal(Form("Alignable entry %s not created. Volume path %s not valid", symname.Data(),volpath.Data()));
+     TGeoPNEntry *alignableEntry = gGeoManager->GetAlignableEntry(symname.Data());
+    const char *path = alignableEntry->GetTitle();
+    if (!gGeoManager->cd(path))
+      AliFatal(Form("Volume path %s not valid!",path));
+    TGeoHMatrix* globMatrix = gGeoManager->GetCurrentMatrix();
+    TGeoHMatrix* matTtoL = fTPCParam->Tracking2LocalMatrix(globMatrix,36+18+cnt-1);
+    alignableEntry->SetMatrix(matTtoL);
     modnum++;
   }
 }