]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - STEER/AliSimulation.cxx
Set beam energy in GRP object to energy per charge.
[u/mrichter/AliRoot.git] / STEER / AliSimulation.cxx
index 0bb13bcbcff28ae85d1f9e48fd437e20c6da8f85..ca5b3a6948c6a5f793618bbb2a89ef1ca7009a8b 100644 (file)
 #include "AliAlignObj.h"
 #include "AliCDBEntry.h"
 #include "AliCDBManager.h"
+#include "AliGRPManager.h"
 #include "AliCDBStorage.h"
 #include "AliCTPRawData.h"
 #include "AliCentralTrigger.h"
@@ -191,6 +192,10 @@ AliSimulation::AliSimulation(const char* configFileName,
   fSetRunNumberFromDataCalled(kFALSE),
   fEmbeddingFlag(kFALSE),
   fLego(NULL),
+  fKey(0),
+  fUseVertexFromCDB(0),
+  fUseMagFieldFromGRP(0),
+  fGRPWriteLocation(Form("local://%s", gSystem->pwd())),
   fQADetectors("ALL"),                  
   fQATasks("ALL"),     
   fRunQA(kTRUE), 
@@ -362,7 +367,8 @@ void AliSimulation::SetCDBLock() {
   // Set CDB lock: from now on it is forbidden to reset the run number
   // or the default storage or to activate any further storage!
   
-  AliCDBManager::Instance()->SetLock(1);
+  ULong_t key = AliCDBManager::Instance()->SetLock(1);
+  if (key) fKey = key;
 }
 
 //_____________________________________________________________________________
@@ -645,9 +651,7 @@ Bool_t AliSimulation::Run(Int_t nEvents)
   if (!AliGeomManager::GetGeometry()) {
     // Initialize the geometry manager
     AliGeomManager::LoadGeometry("geometry.root");
-  AliSysInfo::AddStamp("GetGeometry");
-
-    
+    AliSysInfo::AddStamp("GetGeometry");
 //    // Check that the consistency of symbolic names for the activated subdetectors
 //    // in the geometry loaded by AliGeomManager
 //    AliRunLoader* runLoader = LoadRun("READ");
@@ -983,13 +987,50 @@ Bool_t AliSimulation::RunSimulation(Int_t nEvents)
   //
   // Initialize ALICE Simulation run
   //
-
   gAlice->Announce();
 
+  //
+  // If requested set the mag. field from the GRP entry.
+  // After this the field is loccked and cannot be changed by Config.C
+  if (fUseMagFieldFromGRP) {
+      AliGRPManager grpM;
+      grpM.ReadGRPEntry();
+      grpM.SetMagField();
+      AliInfo("Field is locked now. It cannot be changed in Config.C");
+  }
+//
+// Execute Config.C
   gROOT->LoadMacro(fConfigFileName.Data());
   gInterpreter->ProcessLine(gAlice->GetConfigFunction());
   AliSysInfo::AddStamp("RunSimulation_Config");
 
+//
+// If requested obtain the vertex position and vertex sigma_z from the CDB
+// This overwrites the settings from the Config.C  
+  if (fUseVertexFromCDB) {
+      Double_t vtxPos[3] = {0., 0., 0.}; 
+      Double_t vtxSig[3] = {0., 0., 0.};
+      AliCDBEntry* entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertex");
+      AliESDVertex* vertex = dynamic_cast<AliESDVertex*> (entry->GetObject());
+      if(vertex) {
+         if(vertex->GetXRes()>2.8) { // > pipe radius --> it's a dummy object, don't use it 
+             entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertexSPD");
+             vertex = dynamic_cast<AliESDVertex*> (entry->GetObject());
+         }
+      }
+      if (vertex) {
+         vertex->GetXYZ(vtxPos);
+         vertex->GetSigmaXYZ(vtxSig);
+         AliInfo("Overwriting Config.C vertex settings !");
+         AliInfo(Form("Vertex position from OCDB entry: x = %13.3f, y = %13.3f, z = %13.3f (sigma = %13.3f)\n",
+                vtxPos[0], vtxPos[1], vtxPos[2], vtxSig[2]));
+
+         AliGenerator *gen = gAlice->GetMCApp()->Generator();
+         gen->SetOrigin(vtxPos[0], vtxPos[1], vtxPos[2]);   // vertex position
+         gen->SetSigmaZ(vtxSig[2]);
+      }
+  }
+  
   if(AliCDBManager::Instance()->GetRun() >= 0) { 
     AliRunLoader::Instance()->SetRunNumber(AliCDBManager::Instance()->GetRun());
     AliRunLoader::Instance()->SetNumberOfEventsPerRun(fNEvents);
@@ -998,7 +1039,8 @@ Bool_t AliSimulation::RunSimulation(Int_t nEvents)
   }
   
    AliRunLoader::Instance()->CdGAFile();
-    
+   
+
    AliPDG::AddParticlesToPdgDataBase();  
 
    gMC->SetMagField(TGeoGlobalMagField::Instance()->GetField());
@@ -1101,8 +1143,7 @@ Bool_t AliSimulation::RunSimulation(Int_t nEvents)
     loaderName += "Loader";
     AliLoader* loader = runLoader->GetLoader(loaderName);
     if (!loader) {
-      AliError(Form("RunSimulation", "no loader for %s found\n"
-                    "Number of events per file not set for %s %s", 
+      AliError(Form("RunSimulation no loader for %s found\n Number of events per file not set for %s %s", 
                     detName, typeName, detName));
       continue;
     }
@@ -1341,6 +1382,7 @@ Bool_t AliSimulation::WriteRawData(const char* detectors,
        sprintf(command, "rm -r raw%d", iEvent);
        gSystem->Exec(command);
       }
+      delete runLoader;
     }
 
     if (rootOutput) {
@@ -1750,7 +1792,7 @@ Bool_t AliSimulation::IsSelected(TString detName, TString& detectors) const
 }
 
 //_____________________________________________________________________________
-Bool_t AliSimulation::ConvertRaw2SDigits(const char* rawDirectory, const char* esdFileName) 
+Bool_t AliSimulation::ConvertRaw2SDigits(const char* rawDirectory, const char* esdFileName, Int_t N
 {
 //
 // Steering routine  to convert raw data in directory rawDirectory/ to fake SDigits. 
@@ -1795,6 +1837,7 @@ Bool_t AliSimulation::ConvertRaw2SDigits(const char* rawDirectory, const char* e
    AliRunLoader::Instance()->LoadKinematics("RECREATE");
    AliRunLoader::Instance()->LoadTrackRefs("RECREATE");
    AliRunLoader::Instance()->LoadHits("all","RECREATE");
+
    //
    // Save stuff at the beginning of the file to avoid file corruption
    AliRunLoader::Instance()->CdGAFile();
@@ -1811,11 +1854,17 @@ Bool_t AliSimulation::ConvertRaw2SDigits(const char* rawDirectory, const char* e
     AliRunLoader* runLoader = AliRunLoader::Instance();
     //
     // Open esd file if available
-    TFile* esdFile = TFile::Open(esdFileName);
+    TFile* esdFile = 0;
     TTree* treeESD = 0;
-    AliESDEvent* esd = new AliESDEvent();
-    esdFile->GetObject("esdTree", treeESD);
-    if (treeESD) esd->ReadFromTree(treeESD);
+    AliESDEvent* esd = 0;
+    if (esdFileName && (strlen(esdFileName)>0)) {
+      esdFile = TFile::Open(esdFileName);
+      if (esdFile) {
+        esd = new AliESDEvent();
+        esdFile->GetObject("esdTree", treeESD);
+        if (treeESD) esd->ReadFromTree(treeESD);
+      }
+    }
 
     //
     // Create the RawReader
@@ -1836,24 +1885,29 @@ Bool_t AliSimulation::ConvertRaw2SDigits(const char* rawDirectory, const char* e
     //
     // Get Header
     AliHeader* header = runLoader->GetHeader();
-    //
-    TString detStr = fMakeSDigits;
     // Event loop
     Int_t nev = 0;
     while(kTRUE) {
        if (!(rawReader->NextEvent())) break;
+       runLoader->SetEventNumber(nev);
+        runLoader->GetHeader()->Reset(rawReader->GetRunNumber(), 
+                                      nev, nev);
+       runLoader->TreeE()->Fill();
+        runLoader->GetEvent(nev);
+        AliInfo(Form("We are at event %d",nev));
        //
        // Detector loop
+        TString detStr = fMakeSDigits;
        for (iDet = 0; iDet < detArray->GetEntriesFast(); iDet++) {
            AliModule* det = (AliModule*) detArray->At(iDet);
            if (!det || !det->IsActive()) continue;
            if (IsSelected(det->GetName(), detStr)) {
-             AliInfo(Form("Calling Raw2SDigits for %s\n", det->GetName()));
+             AliInfo(Form("Calling Raw2SDigits for %s", det->GetName()));
              det->Raw2SDigits(rawReader);
              rawReader->Reset();
            }
        } // detectors
-
+        
 
        //
        //  If ESD information available obtain reconstructed vertex and store in header.
@@ -1869,13 +1923,14 @@ Bool_t AliSimulation::ConvertRaw2SDigits(const char* rawDirectory, const char* e
            mcHeader->SetPrimaryVertex(mcV);
            header->Reset(0,nev);
            header->SetGenEventHeader(mcHeader);
-           printf("***** Saved vertex %f %f %f \n", position[0], position[1], position[2]);
+           AliInfo(Form("***** Saved vertex %f %f %f \n", position[0], position[1], position[2]));
        }
-       nev++;
 //
 //      Finish the event
-       runLoader->TreeE()->Fill();
-       runLoader->SetNextEvent();
+        AliInfo(Form("Finished event %d",nev));
+       nev++;
+        if (N>0&&nev>=N)
+          break;
     } // events
  
     delete rawReader;
@@ -2083,7 +2138,7 @@ Bool_t AliSimulation::RunQA()
 {
        // run the QA on summable hits, digits or digits
        
-  if(!gAlice) return kFALSE;
+    //if(!gAlice) return kFALSE;
        AliQAManager::QAManager()->SetRunLoader(AliRunLoader::Instance()) ;
 
        TString detectorsw("") ;  
@@ -2215,10 +2270,11 @@ void AliSimulation::WriteGRPEntry()
   grpObj->SetBeamEnergyIsSqrtSHalfGeV(); // new format of GRP: store sqrt(s)/2 in GeV
 
   const AliGenerator *gen = gAlice->GetMCApp()->Generator();
+  Int_t a = 0;
+  Int_t z = 0;
+
   if (gen) {
-    grpObj->SetBeamEnergy(gen->GetEnergyCMS()/2);
     TString projectile;
-    Int_t a,z;
     gen->GetProjectile(projectile,a,z);
     TString target;
     gen->GetTarget(target,a,z);
@@ -2226,9 +2282,15 @@ void AliSimulation::WriteGRPEntry()
     beamType.ReplaceAll(" ","");
     if (!beamType.CompareTo("-")) {
       grpObj->SetBeamType("UNKNOWN");
+      grpObj->SetBeamEnergy(gen->GetEnergyCMS()/2);
     }
     else {
       grpObj->SetBeamType(beamType);
+      if (z != 0) {
+         grpObj->SetBeamEnergy(gen->GetEnergyCMS()/2 * a / z);
+      } else {
+         grpObj->SetBeamEnergy(gen->GetEnergyCMS()/2 );
+      }
       // Heavy ion run, the event specie is set to kHighMult
       fEventSpecie = AliRecoParam::kHighMult;
       if ((strcmp(beamType,"p-p") == 0) ||
@@ -2292,6 +2354,11 @@ void AliSimulation::WriteGRPEntry()
 
   // Now store the entry in OCDB
   AliCDBManager* man = AliCDBManager::Instance();
+  
+  man->SetLock(0, fKey);
+  
+  AliCDBStorage* sto = man->GetStorage(fGRPWriteLocation.Data());
+  
 
   AliCDBId id("GRP/GRP/Data", man->GetRun(), man->GetRun(), 1, 1);
   AliCDBMetaData *metadata= new AliCDBMetaData();
@@ -2299,7 +2366,8 @@ void AliSimulation::WriteGRPEntry()
   metadata->SetResponsible("alice-off@cern.ch");
   metadata->SetComment("Automatically produced GRP entry for Monte Carlo");
  
-  man->Put(grpObj,id,metadata);
+  sto->Put(grpObj,id,metadata);
+  man->SetLock(1, fKey);
 }