Last fixes before declaring (almost) success
authorcholm <cholm@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 23 Mar 2006 17:10:48 +0000 (17:10 +0000)
committercholm <cholm@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 23 Mar 2006 17:10:48 +0000 (17:10 +0000)
FMD/AliFMD.cxx
FMD/AliFMDDisplay.cxx
FMD/AliFMDDisplay.h
FMD/AliFMDInput.cxx
FMD/AliFMDInput.h
FMD/AliFMDRawReader.cxx
FMD/AliFMDRawReader.h
FMD/libFMDutil.pkg
FMD/scripts/CheckAlign.C [new file with mode: 0644]
FMD/scripts/CheckRaw.C [new file with mode: 0644]
FMD/scripts/MakeAlignment.C

index 4dceb19..7b598c2 100644 (file)
@@ -511,8 +511,8 @@ void
 AliFMD::Init()
 {
   AliDebug(1, "Initialising FMD detector object");
-  AliFMDGeometry*  fmd = AliFMDGeometry::Instance();
-  fmd->InitTransformations();
+  // AliFMDGeometry*  fmd = AliFMDGeometry::Instance();
+  // fmd->InitTransformations();
 }
 
 //____________________________________________________________________
index 9a804a0..a4c4ed8 100644 (file)
@@ -284,6 +284,13 @@ AliFMDDisplay::ProcessDigit(AliFMDDigit* digit)
 
 //____________________________________________________________________
 Bool_t 
+AliFMDDisplay::ProcessRaw(AliFMDDigit* digit)
+{
+  return ProcessDigit(digit);
+}
+
+//____________________________________________________________________
+Bool_t 
 AliFMDDisplay::ProcessRecPoint(AliFMDRecPoint* recpoint)
 {
   if (!recpoint) { AliError("No recpoint");   return kFALSE; }
index 1d2af93..4b87aa4 100644 (file)
@@ -39,6 +39,7 @@ public:
   virtual Bool_t End();
   virtual Bool_t ProcessHit(AliFMDHit* hit, TParticle* p);
   virtual Bool_t ProcessDigit(AliFMDDigit* digit);
+  virtual Bool_t ProcessRaw(AliFMDDigit* digit);
   virtual Bool_t ProcessRecPoint(AliFMDRecPoint* recpoint);
   virtual Int_t  LookupColor(Float_t x, Float_t max)  const;
 protected:
index ce67516..a659264 100644 (file)
 #include "AliRunLoader.h"       // ALIRUNLOADER_H
 #include "AliRun.h"             // ALIRUN_H
 #include "AliStack.h"           // ALISTACK_H
+#include "AliRawReaderFile.h"   // ALIRAWREADERFILE_H
 #include "AliFMD.h"             // ALIFMD_H
 #include "AliFMDHit.h"         // ALIFMDHIT_H
 #include "AliFMDDigit.h"       // ALIFMDDigit_H
 #include "AliFMDRecPoint.h"    // ALIFMDRECPOINT_H
+#include "AliFMDRawReader.h"    // ALIFMDRAWREADER_H
 #include <AliESD.h>
 #include <AliESDFMD.h>
 #include <AliCDBManager.h>
@@ -64,6 +66,7 @@ AliFMDInput::AliFMDInput()
     fRun(0), 
     fStack(0),
     fFMDLoader(0), 
+    fReader(0),
     fFMD(0),
     fMainESD(0),
     fESD(0),
@@ -78,6 +81,7 @@ AliFMDInput::AliFMDInput()
     fArrayD(0),
     fArrayS(0), 
     fArrayR(0), 
+    fArrayA(0), 
     fTreeMask(0), 
     fIsInit(kFALSE)
 {
@@ -96,6 +100,7 @@ AliFMDInput::AliFMDInput(const char* gAliceFile)
     fRun(0), 
     fStack(0),
     fFMDLoader(0), 
+    fReader(0),
     fFMD(0),
     fMainESD(0),
     fESD(0),
@@ -110,6 +115,7 @@ AliFMDInput::AliFMDInput(const char* gAliceFile)
     fArrayD(0),
     fArrayS(0), 
     fArrayR(0), 
+    fArrayA(0), 
     fTreeMask(0), 
     fIsInit(kFALSE)
 {
@@ -185,6 +191,12 @@ AliFMDInput::Init()
     fChainE->SetBranchAddress("ESD", &fMainESD);
   }
     
+  if (TESTBIT(fTreeMask, kRaw)) {
+    AliInfo("Getting FMD raw data digits");
+    fArrayA = new TClonesArray("AliFMDDigit");
+    fReader = new AliRawReaderFile(-1);
+  }
+  
   // Optionally, get the geometry 
   if (TESTBIT(fTreeMask, kGeometry)) {
     TString fname(fRun->GetGeometryFileName());
@@ -276,9 +288,7 @@ AliFMDInput::Begin(Int_t event)
     fTreeR = fFMDLoader->TreeR();
     if (!fArrayR) fArrayR = new TClonesArray("AliFMDRecPoint");
     fTreeR->SetBranchAddress("FMD",  &fArrayR);
-  }
-
-  // Possibly load FMD ESD information 
+  }  // Possibly load FMD ESD information 
   if (TESTBIT(fTreeMask, kESD)) {
     AliInfo("Getting FMD event summary data");
     Int_t read = fChainE->GetEntry(event);
@@ -286,6 +296,14 @@ AliFMDInput::Begin(Int_t event)
     fESD = fMainESD->GetFMDData();
     if (!fESD) return kFALSE;
   }
+  // Possibly load FMD Digit information 
+  if (TESTBIT(fTreeMask, kRaw)) {
+    AliInfo("Getting FMD raw data digits");
+    if (!fReader->NextEvent()) return kFALSE;
+    AliFMDRawReader r(fReader, 0);
+    fArrayA->Clear();
+    r.ReadAdcs(fArrayA);
+  }
   
   return kTRUE;
 }
@@ -309,6 +327,8 @@ AliFMDInput::Event()
     if (!ProcessDigits()) return kFALSE;
   if (TESTBIT(fTreeMask, kSDigits)) 
     if (!ProcessSDigits()) return kFALSE;
+  if (TESTBIT(fTreeMask, kRaw)) 
+    if (!ProcessRawDigits()) return kFALSE;
   if (TESTBIT(fTreeMask, kRecPoints)) 
     if (!ProcessRecPoints()) return kFALSE;
   if (TESTBIT(fTreeMask, kESD))
@@ -393,6 +413,22 @@ AliFMDInput::ProcessSDigits()
 
 //____________________________________________________________________
 Bool_t 
+AliFMDInput::ProcessRawDigits()
+{
+  // Read the digit tree, and pass each digit to the member function
+  // ProcessDigit.
+  Int_t nDigit = fArrayA->GetEntries();
+  if (nDigit <= 0) return kTRUE;
+  for (Int_t j = 0; j < nDigit; j++) {
+    AliFMDDigit* digit = static_cast<AliFMDDigit*>(fArrayA->At(j));
+    if (!digit) continue;
+    if (!ProcessRawDigit(digit)) return kFALSE;
+  }    
+  return kTRUE;
+}
+
+//____________________________________________________________________
+Bool_t 
 AliFMDInput::ProcessRecPoints()
 {
   // Read the reconstrcted points tree, and pass each reconstruction
index 84df8e6..3dec389 100644 (file)
@@ -21,6 +21,7 @@ class AliRunLoader;
 class AliLoader;
 class AliStack;
 class AliRun;
+class AliRawReader;
 class AliFMD;
 class AliFMDHit;
 class AliFMDDigit;
@@ -47,6 +48,7 @@ public:
     kHeader,          // Header information 
     kRecPoints,       // Reconstructed points
     kESD,             // Load ESD's
+    kRaw,             // Read raw data 
     kGeometry         // Not really a tree 
   };
   AliFMDInput();
@@ -67,11 +69,13 @@ public:
   virtual Bool_t ProcessHits();
   virtual Bool_t ProcessDigits();
   virtual Bool_t ProcessSDigits();
+  virtual Bool_t ProcessRawDigits();
   virtual Bool_t ProcessRecPoints();
 
   virtual Bool_t ProcessHit(AliFMDHit*, TParticle*)  { return kTRUE; }
   virtual Bool_t ProcessDigit(AliFMDDigit*)          { return kTRUE; }
   virtual Bool_t ProcessSDigit(AliFMDSDigit*)        { return kTRUE; }
+  virtual Bool_t ProcessRawDigit(AliFMDDigit*)       { return kTRUE; }
   virtual Bool_t ProcessRecPoint(AliFMDRecPoint*)    { return kTRUE; }
   virtual Bool_t ProcessESD(AliESDFMD*)              { return kTRUE; }
   
@@ -81,6 +85,7 @@ protected:
   AliRun*       fRun;        // Run information
   AliStack*     fStack;      // Stack of particles 
   AliLoader*    fFMDLoader;  // Loader of FMD data 
+  AliRawReader* fReader;     // Raw data reader 
   AliFMD*       fFMD;        // FMD object
   AliESD*       fMainESD;    // ESD Object
   AliESDFMD*    fESD;        // FMD ESD data  
@@ -89,12 +94,14 @@ protected:
   TTree*        fTreeD;      // Digit tree 
   TTree*        fTreeS;      // SDigit tree 
   TTree*        fTreeR;      // RecPoint tree
+  TTree*        fTreeA;      // Raw data tree
   TChain*       fChainE;     // Chain of ESD's
   TClonesArray* fArrayE;     // Event info array
   TClonesArray* fArrayH;     // Hit info array
   TClonesArray* fArrayD;     // Digit info array
   TClonesArray* fArrayS;     // SDigit info array
-  TClonesArray* fArrayR;     // Mult (single) info array
+  TClonesArray* fArrayR;     // Rec points info array
+  TClonesArray* fArrayA;     // Raw data (digits) info array
   TGeoManager*  fGeoManager; // Geometry manager
   Int_t         fTreeMask;   // Which tree's to load
   Bool_t        fIsInit;
@@ -133,6 +140,15 @@ public:
 };
 
 //____________________________________________________________________
+class AliFMDInputRaw : public AliFMDInput 
+{
+public:
+  AliFMDInputRaw(const char* file="galice.root") 
+    : AliFMDInput(file) { AddLoad(kRaw); }
+  ClassDef(AliFMDInputRaw, 0);
+};
+
+//____________________________________________________________________
 class AliFMDRecPoint;
 class AliFMDInputRecPoints : public AliFMDInput 
 {
index 96fee2f..9b6bccd 100644 (file)
@@ -75,20 +75,19 @@ AliFMDRawReader::AliFMDRawReader(AliRawReader* reader, TTree* tree)
   // Default CTOR
 }
 
-
 //____________________________________________________________________
-void
-AliFMDRawReader::Exec(Option_t*) 
+Bool_t
+AliFMDRawReader::ReadAdcs(TClonesArray* array) 
 {
   // Read raw data into the digits array, using AliFMDAltroReader. 
+  if (!array) {
+    AliError("No TClonesArray passed");
+    return kFALSE;
+  }
   if (!fReader->ReadHeader()) {
-    Error("ReadAdcs", "Couldn't read header");
-    return;
+    AliError("Couldn't read header");
+    return kFALSE;
   }
-
-  TClonesArray* array = new TClonesArray("AliFMDDigit");
-  fTree->Branch("FMD", &array);
-
   // Get sample rate 
   AliFMDParameters* pars = AliFMDParameters::Instance();
 
@@ -129,7 +128,7 @@ AliFMDRawReader::Exec(Option_t*)
       // Loop over the `timebins', and make the digits
       for (size_t i = 0; i < last; i++) {
        if (i < preSamp) continue;
-       Int_t n = array->GetEntries();
+       Int_t    n      = array->GetEntries();
        UShort_t curStr = str + stripMin + i / rate;
        if ((curStr-str) > stripMax) {
          AliError(Form("Current strip is %d but DB says max is %d", 
@@ -146,7 +145,25 @@ AliFMDRawReader::Exec(Option_t*)
        if (r.IsBof()) break;
     }
   } while (true);
-  AliDebug(1, Form("Got a grand total of %d digits", array->GetEntries()));
+  return kTRUE;
+}
+
+  
+
+//____________________________________________________________________
+void
+AliFMDRawReader::Exec(Option_t*) 
+{
+  TClonesArray* array = new TClonesArray("AliFMDDigit");
+  if (!fTree) {
+    AliError("No tree");
+    return;
+  }
+  fTree->Branch("FMD", &array);
+  ReadAdcs(array);
+  Int_t nWrite = fTree->Fill();
+  AliDebug(1, Form("Got a grand total of %d digits, wrote %d bytes to tree", 
+                  array->GetEntries(), nWrite));
 }
 
 #if 0
index 01d9582..aad6491 100644 (file)
@@ -21,6 +21,7 @@
 //____________________________________________________________________
 class AliRawReader;
 class TTree;
+class TClonesArray;
 
 
 //____________________________________________________________________
@@ -29,7 +30,8 @@ class AliFMDRawReader : public TTask
 public:
   AliFMDRawReader(AliRawReader* reader, TTree* array);
   virtual ~AliFMDRawReader() {}
-  virtual void Exec(Option_t* option="");
+  virtual void   Exec(Option_t* option="");
+  virtual Bool_t ReadAdcs(TClonesArray* array);
 protected:
   TTree*        fTree;       //! Pointer to tree to read into 
   AliRawReader* fReader;     //! Pointer to raw reader 
index 6974165..b604c20 100644 (file)
@@ -8,6 +8,7 @@ SRCS            =  AliFMDInput.cxx              \
                   AliFMDAlignFaker.cxx         
 HDRS           =  $(SRCS:.cxx=.h) 
 DHDR           := FMDutilLinkDef.h
+EINCLUDE       := $(ALICE)/RAW
 
 #
 # EOF
diff --git a/FMD/scripts/CheckAlign.C b/FMD/scripts/CheckAlign.C
new file mode 100644 (file)
index 0000000..87218d4
--- /dev/null
@@ -0,0 +1,173 @@
+//____________________________________________________________________
+//
+// $Id$
+//
+// Script that contains a class to compare the raw data written to the
+// digits it's created from.
+//
+// Use the script `Compile.C' to compile this class using ACLic. 
+//
+#include <AliLog.h>
+#include <AliFMDHit.h>
+#include <AliFMDDigit.h>
+#include <AliFMDInput.h>
+#include <AliFMDUShortMap.h>
+#include <AliFMDParameters.h>
+#include <AliFMDGeometry.h>
+#include <AliFMDRing.h>
+#include <AliFMDDetector.h>
+#include <iostream>
+#include <TH3D.h>
+#include <TStyle.h>
+#include <TArrayF.h>
+#include <TParticle.h>
+#include <TCanvas.h>
+
+class CheckAlign : public AliFMDInput
+{
+public:
+  CheckAlign()
+  {
+    AddLoad(kHits);
+    AddLoad(kDigits);
+    AddLoad(kGeometry);
+    AliFMDGeometry* geom = AliFMDGeometry::Instance();
+    geom->Init();
+    // geom->InitTransformations();
+    Double_t iR  = geom->GetRing('I')->GetHighR()+5;
+    Double_t oR  = geom->GetRing('O')->GetHighR()+5;
+    Double_t z1l = geom->GetDetector(1)->GetInnerZ()-5;
+    Double_t z1h = z1l+10;
+    Double_t z2l = geom->GetDetector(2)->GetOuterZ()-5;
+    Double_t z2h = geom->GetDetector(2)->GetInnerZ()+5;
+    Double_t z3l = geom->GetDetector(3)->GetOuterZ()-5;
+    Double_t z3h = geom->GetDetector(3)->GetInnerZ()+5;
+    
+    f1Hits   = new TH3D("hits1", "FMD1 hits", 
+                       300,-iR,iR,300,-iR,iR,100,z1l,z1h);
+    f1Hits->SetMarkerColor(2);
+    f1Hits->SetMarkerStyle(3);
+      
+    f2Hits   = new TH3D("hits2", "FMD2 hits", 
+                       300,-oR,oR,300,-oR,oR,100,z2l,z2h);
+    f2Hits->SetMarkerColor(2);
+    f2Hits->SetMarkerStyle(3);
+
+    f3Hits   = new TH3D("hits3", "FMD3 hits", 
+                       300,-oR,oR,300,-oR,oR,100,z3l,z3h);
+    f3Hits->SetMarkerColor(2);
+    f3Hits->SetMarkerStyle(3);
+
+    f1Digits = new TH3D("digits1", "FMD1 digits", 
+                       300,-iR,iR,300,-iR,iR,100,z1l,z1h); 
+    f1Digits->SetMarkerColor(4);
+    f1Digits->SetMarkerStyle(2);
+
+    f2Digits = new TH3D("digits2", "FMD2 digits", 
+                       300,-oR,oR,300,-oR,oR,100,z2l,z2h);
+    f2Digits->SetMarkerColor(4);
+    f2Digits->SetMarkerStyle(2);
+
+    f3Digits = new TH3D("digits3", "FMD3 hits", 
+                       300,-oR,oR,300,-oR,oR,100,z3l,z3h);
+    f3Digits->SetMarkerColor(4);
+    f3Digits->SetMarkerStyle(2);
+  }
+  Bool_t Init() 
+  {
+    Bool_t ret = AliFMDInput::Init();
+    AliFMDGeometry* geom = AliFMDGeometry::Instance();
+    geom->Init();
+    geom->InitTransformations();
+    AliFMDParameters* param = AliFMDParameters::Instance();
+    param->Init();
+    return ret;
+  }
+  
+  Bool_t ProcessHit(AliFMDHit* hit, TParticle*)
+  {
+    // Cache the energy loss 
+    if (!hit) return kFALSE;
+    UShort_t det = hit->Detector();
+    UShort_t str = hit->Strip();
+    if (str > 511) {
+      AliWarning(Form("Bad strip number %d in hit", str));
+      return kTRUE;
+    }
+    switch (det) {
+    case 1: f1Hits->Fill(hit->X(), hit->Y(), hit->Z()); break;
+    case 2: f2Hits->Fill(hit->X(), hit->Y(), hit->Z()); break;
+    case 3: f3Hits->Fill(hit->X(), hit->Y(), hit->Z()); break;
+    }
+    return kTRUE;
+  }
+  Bool_t ProcessDigit(AliFMDDigit* digit)
+  {
+    // Cache the energy loss 
+    if (!digit) return kFALSE;
+    UShort_t det = digit->Detector();
+    Char_t   rng = digit->Ring();
+    UShort_t sec = digit->Sector();
+    UShort_t str = digit->Strip();
+    if (str > 511) {
+      AliWarning(Form("Bad strip number %d in digit", str));
+      return kTRUE;
+    }
+    AliFMDParameters* param = AliFMDParameters::Instance();
+    if (digit->Counts() < (param->GetPedestal(det, rng, sec, str) + 4 *
+                          param->GetPedestalWidth(det, rng, sec, str)))
+      return kTRUE;
+    AliFMDGeometry* geom = AliFMDGeometry::Instance();
+    Double_t x, y, z;
+    geom->Detector2XYZ(det, rng, sec, str, x, y, z);
+    switch (det) {
+    case 1: f1Digits->Fill(x, y , z); break;
+    case 2: f2Digits->Fill(x, y , z); break;
+    case 3: f3Digits->Fill(x, y , z); break;
+    }
+    return kTRUE;
+  }
+  //__________________________________________________________________
+  Bool_t Finish()
+  {
+    gStyle->SetPalette(1);
+    gStyle->SetOptTitle(0);
+    gStyle->SetCanvasColor(0);
+    gStyle->SetCanvasBorderSize(0);
+    gStyle->SetPadColor(0);
+    gStyle->SetPadBorderSize(0);
+
+    TCanvas* c1 = new TCanvas("FMD1","FMD1");
+    c1->cd();
+    f1Hits->Draw();
+    f1Digits->Draw("same");
+
+    TCanvas* c2 = new TCanvas("FMD2","FMD2");
+    c2->cd();
+    f2Hits->Draw();
+    f2Digits->Draw("same");
+
+    TCanvas* c3 = new TCanvas("FMD3","FMD3");
+    c3->cd();
+    f3Hits->Draw();
+    f3Digits->Draw("same");
+    return kTRUE;
+  }
+protected:
+  TH3D* f1Hits;
+  TH3D* f2Hits;
+  TH3D* f3Hits;
+  TH3D* f1Digits;
+  TH3D* f2Digits;
+  TH3D* f3Digits;
+};
+
+
+//____________________________________________________________________
+//
+// EOF
+//
+
+
+  
+  
diff --git a/FMD/scripts/CheckRaw.C b/FMD/scripts/CheckRaw.C
new file mode 100644 (file)
index 0000000..745a9f5
--- /dev/null
@@ -0,0 +1,88 @@
+//____________________________________________________________________
+//
+// $Id$
+//
+// Script that contains a class to compare the raw data written to the
+// digits it's created from.
+//
+// Use the script `Compile.C' to compile this class using ACLic. 
+//
+#include <AliLog.h>
+#include <AliFMDDigit.h>
+#include <AliFMDInput.h>
+#include <AliFMDUShortMap.h>
+#include <AliFMDParameters.h>
+#include <iostream>
+#include <TStyle.h>
+#include <TArrayF.h>
+#include <TCanvas.h>
+
+class CheckRaw : public AliFMDInput
+{
+public:
+  CheckRaw()
+  {
+    AddLoad(kDigits);
+    AddLoad(kRaw);
+  }
+  Bool_t Init() 
+  {
+    Bool_t ret = AliFMDInput::Init();
+    // AliFMDGeometry* geom = AliFMDGeometry::Instance();
+    // geom->Init();
+    // geom->InitTransformations();
+    AliFMDParameters* param = AliFMDParameters::Instance();
+    param->Init();
+    return ret;
+  }
+  Bool_t ProcessDigit(AliFMDDigit* digit)
+  {
+    // Cache the energy loss 
+    if (!digit) return kFALSE;
+    UShort_t det = digit->Detector();
+    Char_t   rng = digit->Ring();
+    UShort_t sec = digit->Sector();
+    UShort_t str = digit->Strip();
+    if (str > 511) {
+      AliWarning(Form("Bad strip number %d in digit", str));
+      return kTRUE;
+    }
+    fMap(det, rng, sec, str) = digit->Counts();
+    return kTRUE;
+  }
+  Bool_t ProcessRawDigit(AliFMDDigit* digit)
+  {
+    // Cache the energy loss 
+    if (!digit) return kFALSE;
+    UShort_t det = digit->Detector();
+    Char_t   rng = digit->Ring();
+    UShort_t sec = digit->Sector();
+    UShort_t str = digit->Strip();
+    if (str > 511) {
+      AliWarning(Form("Bad strip number %d in digit", str));
+      return kTRUE;
+    }
+    if (digit->Counts() != fMap(det, rng, sec, str) && 
+       fMap(det, rng, sec, str) != 1024) {
+      AliWarning(Form("Mismatch in digit FMD%d%c[%2d,%3d] %d != %d", 
+                     det, rng, sec, str, digit->Counts(), 
+                     fMap(det, rng, sec, str)));
+      return kTRUE;
+    }
+    AliDebug(1, Form("Raw digit FMD%d%c[%2d,%3D] is good", 
+                    det, rng, sec, str));
+    return kTRUE;
+  }
+protected:
+  AliFMDUShortMap fMap;
+};
+
+
+//____________________________________________________________________
+//
+// EOF
+//
+
+
+  
+  
index 69f669c..20aaecd 100644 (file)
@@ -15,8 +15,8 @@ MakeAlignment()
   gSystem->Load("libFMDutil.so");
   AliFMDAlignFaker f;
   f.RemoveAlign(AliFMDAlignFaker::kHalves);
-  f.SetSensorDisplacement(0,0,0,0,0,0);
-  f.SetSensorRotation(0,0,0,0,0,0);
+  f.SetSensorDisplacement(0, 0, 0, 0, 0, 0);
+  f.SetSensorRotation(0, 0, 0, 3, 3, 3);
   f.Exec();
 }
 //____________________________________________________________________