Reconstruction of raw data (T.Kuhr)
[u/mrichter/AliRoot.git] / TPC / AliTPCReconstructor.cxx
index 4df397b01b76b9567f7518e67cccfa231499366a..a491a174359bb40e9549e515aa8f1c787cd872de 100644 (file)
@@ -25,6 +25,7 @@
 #include "AliTPCReconstructor.h"
 #include "AliRunLoader.h"
 #include "AliRun.h"
+#include "AliRawReader.h"
 #include "AliTPCclustererMI.h"
 #include "AliTPCtrackerMI.h"
 #include "AliTPCpidESD.h"
@@ -32,6 +33,7 @@
 
 ClassImp(AliTPCReconstructor)
 
+Double_t AliTPCReconstructor::fgCtgRange = 1.05;
 
 //_____________________________________________________________________________
 void AliTPCReconstructor::Reconstruct(AliRunLoader* runLoader) const
@@ -76,6 +78,42 @@ void AliTPCReconstructor::Reconstruct(AliRunLoader* runLoader) const
   loader->UnloadDigits();
 }
 
+//_____________________________________________________________________________
+void AliTPCReconstructor::Reconstruct(AliRunLoader* runLoader,
+                                     AliRawReader* rawReader) const
+{
+// reconstruct clusters from raw data
+
+  AliLoader* loader = runLoader->GetLoader("TPCLoader");
+  if (!loader) {
+    Error("Reconstruct", "TPC loader not found");
+    return;
+  }
+  loader->LoadRecPoints("recreate");
+
+  AliTPCParam* param = GetTPCParam(runLoader);
+  if (!param) return;
+  AliTPCclustererMI clusterer(param);
+
+  Int_t iEvent = 0;
+  while (rawReader->NextEvent()) {
+    runLoader->GetEvent(iEvent++);
+
+    TTree* treeClusters = loader->TreeR();
+    if (!treeClusters) {
+      loader->MakeTree("R");
+      treeClusters = loader->TreeR();
+    }
+
+    clusterer.SetOutput(treeClusters);
+    clusterer.Digits2Clusters(rawReader);
+         
+    loader->WriteRecPoints("OVERWRITE");
+  }
+
+  loader->UnloadRecPoints();
+}
+
 //_____________________________________________________________________________
 AliTracker* AliTPCReconstructor::CreateTracker(AliRunLoader* runLoader) const
 {