Added possibility to monitor performance in simulation with respect to geometry and...
authoragheata <agheata@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 30 Jan 2013 14:18:02 +0000 (14:18 +0000)
committeragheata <agheata@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 30 Jan 2013 14:18:02 +0000 (14:18 +0000)
STEER/CMakelibSTEER.pkg
STEER/STEER/AliMC.cxx
STEER/STEER/AliMC.h
STEER/STEER/AliSimulation.cxx
STEER/STEER/AliSimulation.h
STEER/STEER/AliTransportMonitor.cxx [new file with mode: 0644]
STEER/STEER/AliTransportMonitor.h [new file with mode: 0644]
STEER/STEERLinkDef.h

index f17f570..6cfc349 100644 (file)
@@ -115,6 +115,7 @@ set ( SRCS
     STEER/AliLHCClockPhase.cxx 
     STEER/AliLTUConfig.cxx
     STEER/AliTriggerUtils.cxx
+    STEER/AliTransportMonitor.cxx
     )
 
 string(REPLACE ".cxx" ".h" HDRS  "${SRCS}")
index 83283e8..1a988ea 100644 (file)
@@ -53,6 +53,7 @@
 #include "AliSimulation.h"
 #include "AliStack.h"
 #include "AliTrackReference.h"
+#include "AliTransportMonitor.h"
 
 using std::endl;
 using std::cout;
@@ -64,6 +65,7 @@ AliMC::AliMC() :
   fSaveRndmStatus(kFALSE),
   fSaveRndmEventStatus(kFALSE),
   fReadRndmStatus(kFALSE),
+  fUseMonitoring(kFALSE),
   fRndmFileName("random.root"),
   fEventEnergy(0),
   fSummEnergy(0),
@@ -75,6 +77,7 @@ AliMC::AliMC() :
   fDecayPdg(0),
   fImedia(0),
   fTransParName("\0"),
+  fMonitor(0),
   fHitLists(0),
   fTmpTreeTR(0),
   fTmpFileTR(0),
@@ -93,6 +96,7 @@ AliMC::AliMC(const char *name, const char *title) :
   fSaveRndmStatus(kFALSE),
   fSaveRndmEventStatus(kFALSE),
   fReadRndmStatus(kFALSE),
+  fUseMonitoring(kFALSE),
   fRndmFileName("random.root"),
   fEventEnergy(0),
   fSummEnergy(0),
@@ -104,6 +108,7 @@ AliMC::AliMC(const char *name, const char *title) :
   fDecayPdg(0),
   fImedia(new TArrayI(1000)),
   fTransParName("\0"),
+  fMonitor(0),
   fHitLists(new TList()),
   fTmpTreeTR(0),
   fTmpFileTR(0),
@@ -125,6 +130,7 @@ AliMC::~AliMC()
   delete fGenerator;
   delete fImedia;
   delete fHitLists;
+  delete fMonitor;
   // Delete track references
 }
 
@@ -456,6 +462,12 @@ void AliMC::FinishRun()
   // Clean generator information
   AliDebug(1, "fGenerator->FinishRun()");
   fGenerator->FinishRun();
+  
+  // Monitoring information
+  if (fMonitor) {
+    fMonitor->Print();
+    fMonitor->Export("timing.root");
+  }  
 
   //Output energy summary tables
   AliDebug(1, "EnergySummary()");
@@ -510,8 +522,25 @@ void AliMC::Stepping()
       FixParticleDecaytime();
   } 
     
-
-  
+  // --- If monitoring timing was requested, monitor the step
+  if (fUseMonitoring) {
+    if (!fMonitor) {
+      fMonitor = new AliTransportMonitor(gMC->NofVolumes()+1);
+      fMonitor->Start();
+    }  
+    if (gMC->IsNewTrack() || gMC->TrackTime() == 0. || gMC->TrackStep()<1.1E-10) {
+      fMonitor->DummyStep();
+    } else {
+    // Normal stepping
+      Int_t copy;
+      Int_t volId = gMC->CurrentVolID(copy);
+      Int_t pdg = gMC->TrackPid();
+      TLorentzVector xyz, pxpypz;
+      gMC->TrackPosition(xyz);
+      gMC->TrackMomentum(pxpypz);
+      fMonitor->StepInfo(volId, pdg, pxpypz.E(), xyz.X(), xyz.Y(), xyz.Z());
+    }  
+  }
   //
   // --- If lego option, do it and leave 
   if (AliSimulation::Instance()->Lego())
index ed1a391..b35cd8f 100644 (file)
@@ -24,6 +24,7 @@ class TTree;
 
 class AliGenerator;
 class AliTrackReference;
+class AliTransportMonitor;
 
 class AliMC : public TVirtualMCApplication {
 public:
@@ -110,7 +111,9 @@ public:
 // Geometry related
    void           SetGeometryFromCDB();
    Bool_t         IsGeometryFromCDB() const;
-   
+// Monitor transport   
+   void           SetUseMonitoring(Bool_t flag=kTRUE)      { fUseMonitoring = flag; }
+   AliTransportMonitor *GetTransportMonitor() const        { return fMonitor; }
 // Random number generator status
    void           SetSaveRndmStatus(Bool_t value)          { fSaveRndmStatus = value; }  
    void           SetSaveRndmStatusPerEvent(Bool_t value)  { fSaveRndmEventStatus = value; }  
@@ -128,6 +131,7 @@ public:
    Bool_t         fSaveRndmStatus;    //! Options to save random engine status
    Bool_t         fSaveRndmEventStatus; //! Options to save random engine status for each event
    Bool_t         fReadRndmStatus;    //! Options to read random engine status
+   Bool_t         fUseMonitoring;     //! Activate monitoring
    TString        fRndmFileName;      //! The file name of random engine status to be read in
    TArrayF        fEventEnergy;       //! Energy deposit for current event
    TArrayF        fSummEnergy;        //! Energy per event in each volume
@@ -139,6 +143,7 @@ public:
    Int_t          fDecayPdg;          //  PDG code of particle with forced decay length
    TArrayI       *fImedia;            //! Array of correspondence between media and detectors
    TString        fTransParName;      //  Name of the transport parameters file
+   AliTransportMonitor *fMonitor;     //! Transport monitoring tool
    TList         *fHitLists;          //! Lists of hits to be remapped by PurifyKine
    //Temporary Track Reference tree related
    TTree         *fTmpTreeTR;            //! Temporary track reference tree
index 1c41528..c51d2f4 100644 (file)
@@ -182,6 +182,7 @@ AliSimulation::AliSimulation(const char* configFileName,
   fDeleteIntermediateFiles(kFALSE),
   fWriteSelRawData(kFALSE),
   fStopOnError(kFALSE),
+  fUseMonitoring(kFALSE),
   fNEvents(1),
   fConfigFileName(configFileName),
   fGAliceFileName("galice.root"),
@@ -1003,6 +1004,9 @@ Bool_t AliSimulation::RunSimulation(Int_t nEvents)
     AliError("gAlice was already run. Restart aliroot and try again.");
     return kFALSE;
   }
+  
+  // Setup monitoring if requested
+  gAlice->GetMCApp()->SetUseMonitoring(fUseMonitoring);
 
   AliInfo(Form("initializing gAlice with config file %s",
           fConfigFileName.Data()));
index 7137e49..eb005f5 100644 (file)
@@ -70,6 +70,7 @@ public:
   void           SetAlignObjArray(TObjArray *array)
                    {fAlignObjArray = array;
                   fLoadAlignFromCDB = kFALSE;}
+  void           SetUseMonitoring(Bool_t flag=kTRUE) {fUseMonitoring = flag;}
 
   Bool_t         MisalignGeometry(AliRunLoader *runLoader = NULL);
 
@@ -168,6 +169,7 @@ private:
   Bool_t         fDeleteIntermediateFiles; // delete intermediate raw data files
   Bool_t         fWriteSelRawData;    // write detectors raw data in a separate file accoring to the trigger cluster
   Bool_t         fStopOnError;        // stop or continue on errors
+  Bool_t         fUseMonitoring;      // monitor simulation timing per volume
 
   Int_t          fNEvents;            // number of events
   TString        fConfigFileName;     // name of the config file
diff --git a/STEER/STEER/AliTransportMonitor.cxx b/STEER/STEER/AliTransportMonitor.cxx
new file mode 100644 (file)
index 0000000..5068597
--- /dev/null
@@ -0,0 +1,305 @@
+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ *                                                                        *
+ * Author: The ALICE Off-line Project.                                    *
+ * Contributors are mentioned in the code where appropriate.              *
+ *                                                                        *
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  *
+ **************************************************************************/
+
+// Class that can be plugged in the simulation to monitor transport timing per 
+// particle for each geometry volume.
+//
+//  andrei.gheata@cern.ch 
+
+#include "AliTransportMonitor.h"
+
+#include "TDatabasePDG.h"
+#include "TCanvas.h"
+#include "TMath.h"
+#include "TFile.h"
+#include "TROOT.h"
+#include "THashList.h"
+#include "TH2F.h"
+#include "TGeoManager.h"
+#include "AliPDG.h"
+
+ClassImp(AliTransportMonitor)
+ClassImp(AliTransportMonitor::AliTransportMonitorVol)
+ClassImp(AliTransportMonitor::AliTransportMonitorVol::AliPMonData)
+
+typedef AliTransportMonitor::AliTransportMonitorVol::AliPMonData PMonData;
+
+//______________________________________________________________________________
+AliTransportMonitor::AliTransportMonitorVol::AliTransportMonitorVol()
+                    :TNamed(),
+                     fNtypes(0),
+                     fTotalTime(0),
+                     fPData(0),
+                     fTimeRZ(0),
+                     fParticles()
+{
+// Default constructor
+}
+
+//______________________________________________________________________________
+AliTransportMonitor::AliTransportMonitorVol::~AliTransportMonitorVol()
+{
+// Destructor
+  delete [] fPData;
+  delete fTimeRZ;
+}
+
+//______________________________________________________________________________
+void AliTransportMonitor::AliTransportMonitorVol::StepInfo(
+                                    Int_t pdg,
+                                    Double_t energy, 
+                                    Double_t dt,
+                                    Double_t x, Double_t y, Double_t z)
+{
+// This method is called at each N steps to store timing info.
+  PMonData &data = GetPMonData(pdg);
+  data.fEdt += energy*dt;
+  data.fTime += dt;
+  fTotalTime += dt;
+  if (!fTimeRZ) {
+    Bool_t status = TH1::AddDirectoryStatus();
+    TH1::AddDirectory(kFALSE);
+    fTimeRZ = new TH2F("h2rz", "", 100, -5000., 5000, 100, 0., 1000.);
+    TH1::AddDirectory(status);
+  }
+  Double_t r = TMath::Sqrt(x*x+y*y);
+  fTimeRZ->Fill(z,r,dt);
+}
+
+//______________________________________________________________________________
+PMonData &AliTransportMonitor::AliTransportMonitorVol::GetPMonData(Int_t pdg)
+{
+// Retrieve stored monitoring object for a given pdg type. If not existing 
+// create one.
+
+  // The object could have been retrieved from file, in which case we have to 
+  // build the map.
+  PMonData *data;
+  if (fNtypes) {
+    if (fParticles.empty()) {
+      for (Int_t i=0; i<fNtypes; i++) {
+        data = &fPData[i];
+        fParticles[i] = data->fPDG;
+      }
+    }
+    ParticleMapIt_t it = fParticles.find(pdg);
+    if (it == fParticles.end()) {
+      data = &fPData[fNtypes]; 
+      data->fPDG = pdg;
+      fParticles[pdg] = fNtypes++;
+    } else {
+      data = &fPData[it->second];
+    }  
+  } else {
+     TDatabasePDG *pdgDB = TDatabasePDG::Instance();
+     if (!pdgDB->ParticleList()) AliPDG::AddParticlesToPdgDataBase();
+     Int_t size = pdgDB->ParticleList()->GetSize();
+     fPData = new PMonData[size];
+     data = &fPData[fNtypes];
+     data->fPDG = pdg;
+     fParticles[pdg] = fNtypes++;
+  }
+  return *data;
+}
+
+ClassImp(AliTransportMonitor)
+
+//______________________________________________________________________________
+AliTransportMonitor::AliTransportMonitor()
+                    :TObject(),
+                     fTotalTime(0),
+                     fTimer(),
+                     fVolumeMon(0)
+{
+// Default constructor
+}
+
+//______________________________________________________________________________
+AliTransportMonitor::AliTransportMonitor(Int_t nvolumes)
+                    :TObject(),
+                     fTotalTime(0),
+                     fTimer(),
+                     fVolumeMon(0)
+{
+// Default constructor
+  fVolumeMon = new TObjArray(nvolumes);
+  fVolumeMon->SetOwner();
+  for (Int_t i=0; i<nvolumes; i++) {
+    AliTransportMonitorVol *volMon = new AliTransportMonitorVol();
+    if (gGeoManager) volMon->SetName(gGeoManager->GetListOfUVolumes()->At(i)->GetName());
+    fVolumeMon->Add(volMon);
+  }   
+}
+
+//______________________________________________________________________________
+AliTransportMonitor::~AliTransportMonitor()
+{
+// Destructor
+  delete fVolumeMon;
+}
+
+//______________________________________________________________________________
+void AliTransportMonitor::Print(Option_t *volName) const
+{
+// Inspect the timing statistics for a single volume or for all the setup
+  Int_t uid = -1;
+  Int_t ntotal = 0;
+  if (!fVolumeMon || !(ntotal=fVolumeMon->GetEntriesFast())) {
+     Info("Inspect", "Transport monitor is empty !");
+     return;
+  }   
+  if (strlen(volName)) {
+    TString svname = volName;
+    Int_t i = 0;
+    for (i=1; i<ntotal; i++) if (svname == fVolumeMon->At(i)->GetName()) break;
+    if (i==ntotal) {
+       Error("Inspect", "No monitoring info stored for volume %s", volName);
+       return;
+    }
+    uid = i;
+    AliTransportMonitorVol *volMon = (AliTransportMonitorVol*)fVolumeMon->At(uid);
+    Int_t ntypes = volMon->GetNtypes();
+    if (!ntypes) {
+      Info("Inspect", "No particles crossed volume %s", volName);
+      return;
+    }  
+    Double_t *timeperpart = new Double_t[ntypes];
+    Int_t *isort = new Int_t[ntypes];
+    Double_t timepervol = 0.;
+    for (i=0; i<ntypes; i++) {
+      timeperpart[i] = volMon->GetTime(i); 
+      timepervol += timeperpart[i];
+    }
+    printf("Volume %s: Transport time: %g%% of %g [s]\n", volMon->GetName(), 100.*timepervol/fTotalTime, fTotalTime);
+    TMath::Sort(ntypes, timeperpart, isort, kTRUE);
+    TString particle;
+    TDatabasePDG *pdgDB = TDatabasePDG::Instance();    
+    if (!pdgDB->ParticleList()) AliPDG::AddParticlesToPdgDataBase();
+    for (i=0; i<ntypes; i++)  {
+       timeperpart[i] /=  timepervol;
+       particle = pdgDB->GetParticle(volMon->GetPDG(isort[i]))->GetName();
+       printf("   %s: %g%%  mean energy: %g\n", particle.Data(), 100.*timeperpart[i], volMon->GetEmed(i));
+    }
+    if (volMon->GetHistogram()) {
+      TCanvas *c1 = (TCanvas*)gROOT->GetListOfCanvases()->FindObject("crz");
+      if (!c1) c1 = new TCanvas("crz");
+      c1->cd();
+      volMon->GetHistogram()->GetXaxis()->SetTitle("z [cm]");
+      volMon->GetHistogram()->GetYaxis()->SetTitle("r [cm]");
+      volMon->GetHistogram()->SetTitle(Form("RZ plot weighted by time spent in %s",volMon->GetName()));
+      volMon->GetHistogram()->Draw();
+    }  
+    return;
+  }
+  // General view
+  TIter next(fVolumeMon);
+  AliTransportMonitorVol *volMon;
+  Int_t ncrossed = 0;
+  TH1F *hnames = new TH1F("volume timing", "relative volume timing", 3,0,3);
+  hnames->SetStats(0);
+  hnames->SetFillColor(38);
+  hnames->SetBit(TH1::kCanRebin);
+  while ((volMon=(AliTransportMonitorVol*)next())) {
+    if (volMon->GetNtypes()) {
+      hnames->Fill(volMon->GetName(), volMon->GetTotalTime());
+      ncrossed++;
+    }
+  }
+  
+  hnames->LabelsDeflate();
+  hnames->GetXaxis()->LabelsOption(">");
+  
+  TCanvas *c = (TCanvas*)gROOT->GetListOfCanvases()->FindObject("cvol_timing");
+  if (!c) c = new TCanvas("cvol_timing");
+  c->cd();
+  c->SetLogy();
+  hnames->Draw();
+  
+  printf("=============================================================================\n");
+  printf("Effective transport time:  %6.2f minutes\n", fTotalTime/60.);
+  printf("Number of crossed volumes: %d from %d\n", ncrossed, fVolumeMon->GetEntriesFast());
+  printf("=============================================================================\n");  
+}     
+
+//______________________________________________________________________________
+void AliTransportMonitor::DummyStep()
+{
+// Reset timer for zero-length steps
+   fTimer.Stop();
+   fTimer.Start(kTRUE);
+}   
+
+//______________________________________________________________________________
+void AliTransportMonitor::StepInfo( Int_t volId,
+                                    Int_t pdg,
+                                    Double_t energy, 
+                                    Double_t x, Double_t y, Double_t z)
+{
+// This method is called at each N steps to store timing info.
+  fTimer.Stop();
+  Double_t dt = fTimer.RealTime();
+  fTotalTime += dt;
+  AliTransportMonitorVol *volMon = (AliTransportMonitorVol*)fVolumeMon->At(volId);
+  volMon->StepInfo(pdg,energy,dt,x,y,z);
+  fTimer.Start(kTRUE);
+}
+
+//______________________________________________________________________________
+void AliTransportMonitor::Start()
+{
+// Start collecting timing information
+  if (fTotalTime > 0) {
+    Info("Start", "Cannot start twice");
+    return;
+  }
+  fTimer.Start(kTRUE);
+}
+
+//______________________________________________________________________________
+void AliTransportMonitor::Stop()
+{
+// Stop the global timer
+  fTimer.Stop();
+}
+
+//______________________________________________________________________________
+void AliTransportMonitor::Export(const char *fname)
+{
+// Export information to file
+  TFile *file = TFile::Open(fname, "RECREATE");
+  Write();
+  file->Write();
+  file->Close();
+}  
+
+//______________________________________________________________________________
+AliTransportMonitor *AliTransportMonitor::Import(const char *fname)
+{
+// Import information from a file
+  TFile *file = TFile::Open(fname);
+  if (!file) {
+    ::Error("Import", "File %s could not be opened", fname);
+    return 0;
+  }
+  AliTransportMonitor *mon = (AliTransportMonitor *)file->Get("AliTransportMonitor");
+  if (!mon) {
+    ::Error("Import", "No AliTransportMonitor object found n file %s", fname);
+    return 0;
+  }
+  AliPDG::AddParticlesToPdgDataBase();
+  return mon;
+}
+
+    
diff --git a/STEER/STEER/AliTransportMonitor.h b/STEER/STEER/AliTransportMonitor.h
new file mode 100644 (file)
index 0000000..7300016
--- /dev/null
@@ -0,0 +1,98 @@
+#ifndef ALI_TRANSPORTMONITOR__H
+#define ALI_TRANSPORTMONITOR__H
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice                               */
+
+#ifndef ROOT_TNamed
+#include "TNamed.h"
+#endif
+
+#include <map>
+
+#ifndef ROOT_TStopwatch
+#include "TStopwatch.h"
+#endif
+
+// Class that can be plugged in the simulation to monitor transport timing per 
+// particle for each geometry volume.
+//
+//  andrei.gheata@cern.ch 
+
+class TH2F;
+
+
+//______________________________________________________________________________
+class AliTransportMonitor : public TObject {
+public:
+  //________________________________________________________________
+  class AliTransportMonitorVol : public TNamed {
+  public:
+    //___________________________________________________
+    class AliPMonData {
+    public:
+      Int_t         fPDG;        // particle PDG
+      Double_t      fEdt;        // Energy * dt integral
+      Double_t      fTime;       // Total transport time for the particle in this volume
+      AliPMonData() : fPDG(0), fEdt(0), fTime(0) {}
+      virtual ~AliPMonData() {}
+      ClassDef(AliPMonData, 1)     // Basic monitoring info structure
+    };   
+    //___________________________________________________   
+    AliTransportMonitorVol();
+    virtual ~AliTransportMonitorVol();
+    
+    Int_t           GetNtypes() const {return fNtypes;}
+    void            StepInfo(Int_t    pdg,
+                             Double_t energy, 
+                             Double_t dt,
+                             Double_t x, Double_t y, Double_t z);
+    Double_t        GetTotalTime() const       {return fTotalTime;}
+    Double_t        GetTime(Int_t itype) const {return fPData[itype].fTime;}
+    Double_t        GetEmed(Int_t itype) const {return (fTotalTime>0)?fPData[itype].fEdt/fTotalTime : 0.;}
+    Int_t           GetPDG(Int_t itype)  const {return fPData[itype].fPDG;}
+    TH2F           *GetHistogram() const {return fTimeRZ;}
+    private:
+      AliPMonData  &GetPMonData(Int_t pdg);
+      AliTransportMonitorVol(const AliTransportMonitorVol& other) : TNamed(other), fNtypes(0), fTotalTime(0), fPData(0), fTimeRZ(0), fParticles() {}
+      AliTransportMonitorVol &operator=(const AliTransportMonitorVol&) {return *this;}
+    private:
+      Int_t         fNtypes;     // Number of different particle types
+      Double_t      fTotalTime;  // Total time spent in this volume
+      AliPMonData  *fPData;      //[fNtypes] Array of particle data
+      TH2F         *fTimeRZ;     // Timing R-Z histogram per volume
+      typedef std::map<Int_t, Int_t>         ParticleMap_t;
+      typedef ParticleMap_t::iterator        ParticleMapIt_t;
+      ParticleMap_t fParticles;  //! Map of stored particla data  
+  
+  ClassDef(AliTransportMonitorVol,1)  // Helper to hold particle info per volume
+  };
+  //________________________________________________________________
+private:
+  AliTransportMonitor(const AliTransportMonitor&other) : TObject(other), fTotalTime(0), fTimer(), fVolumeMon(0) {}
+  AliTransportMonitor &operator=(const AliTransportMonitor&) {return *this;}
+public:
+  AliTransportMonitor();
+  AliTransportMonitor(Int_t nvolumes);
+  virtual ~AliTransportMonitor();
+  
+  void              StepInfo(Int_t    volId,
+                             Int_t    pdg,
+                             Double_t energy, 
+                             Double_t x, Double_t y, Double_t z);
+  void              Print(Option_t *volName="") const;
+  void              DummyStep();
+  void              Start();
+  void              Stop();
+  void              Export(const char *fname);
+  static AliTransportMonitor *Import(const char *fname);
+private:
+  Double_t          fTotalTime;  // Total simulation time
+  TStopwatch        fTimer;      //! Global timer
+  TObjArray        *fVolumeMon;  // Array of monitoring objects per volume
+  
+ClassDef(AliTransportMonitor,1)  // Class to monitor timing per volume 
+};
+//______________________________________________________________________________
+
+
+#endif //ALI_TRANSPORTMONITOR__H
index 5632a44..5a8bc6e 100644 (file)
@@ -7,7 +7,7 @@
 #pragma link off all globals;
 #pragma link off all classes;
 #pragma link off all functions;
+
 #pragma link C++ global gAlice;
 #pragma link C++ global gMC;
  
 
 #pragma link C++ class AliLTUConfig+;
 
+#pragma link C++ class AliTransportMonitor+;
+#pragma link C++ class AliTransportMonitor::AliTransportMonitorVol+;
+#pragma link C++ struct AliTransportMonitor::AliTransportMonitorVol::AliPMonData+;
+
 #pragma link C++ typedef AliLHCDipValD;         
 #pragma link C++ typedef AliLHCDipValI;         
 #pragma link C++ typedef AliLHCDipValF;         
 #pragma link C++ typedef AliLHCDipValC;
+
 #endif