emc_complete_corr.C File Reference

Go to the source code of this file.

Functions
	emc_complete_corr (Int_t nEvents=10, Float_t mom=1.)

Function Documentation

emc_complete_corr	(	Int_t	nEvents = 10,
		Float_t	mom = 1.
	)

Modified 17.12.08 by M. Al-Turany The hit producing is added to the simulation session, so that the EmcHits are produced on the fly during the simulation

Set the digitization parameters

Get the run time data base for this session and set the needed input

Parameters created for this simulation goes to the out put

Add Hit producer task to the simulation

Initialize the session

After initialization now we can save the field parameters

All parameters are initialized and ready to be saved

Definition at line 9 of file emc_complete_corr.C.

References Bool_t, boxGen, Cave, ctime, digiFile, Double_t, Emc, emcDigiFile, emcHdrFiller, emcHitProd, emcHitsToWaveform, emcMakeBump, emcMakeCluster, emcMakeRecoHit, emcWaveformToDigi, fField, fRun, gDebug, kParameterMerged, Magnet, nEvents, output, Par, parFile, parIo1, Pipe, primGen, printf(), rootlogon(), rtdb, rtime, PndEmc::SetGeometryVersion(), PndMultiFieldPar::SetParameters(), PndEmcMakeCluster::SetStorageOfData(), PndEmcMakeBump::SetStorageOfData(), PndEmcHitsToWaveform::SetStorageOfData(), PndEmcWaveformToDigi::SetStorageOfData(), PndEmcMakeCorr::SetStorageOfData(), PndEmcHitProducer::SetStorageOfData(), PndEmc::SetStorageOfData(), timer, and TString.

                                                       {
  TStopwatch timer;
  timer.Start();
  gDebug=0;
  // Load basic libraries
  // If it does not work,  please check the path of the libs and put it by hands
  gROOT->LoadMacro("$VMCWORKDIR/gconfig/rootlogon.C");
  rootlogon();
                
  FairRunSim *fRun = new FairRunSim();
        
  TString digiFile = "emc.par";
        
  TString parFile = "simparams.root";
        
  // set the MC version used
  // ------------------------
        
   //fRun->SetName("TGeant3");
   fRun->SetName("TGeant4");
        
  TString transportModel;
  
  transportModel = fRun->GetName();
  cout<<"transport model = "<<transportModel <<endl;
  
  fRun->SetOutputFile("emc_complete_corr.root");


  TString emcDigiFile = gSystem->Getenv("VMCWORKDIR");
  emcDigiFile += "/macro/params/";
  emcDigiFile += digiFile;
  FairRuntimeDb* rtdb = fRun->GetRuntimeDb();
  FairParAsciiFileIo* parIo1 = new FairParAsciiFileIo();
  parIo1->open(emcDigiFile.Data(),"in");
  rtdb->setFirstInput(parIo1);        
  Bool_t kParameterMerged=kTRUE;

  FairParRootFileIo* output=new FairParRootFileIo(kParameterMerged);
  output->open("simparams.root");
  rtdb->setOutput(output);
        
        
  // Set Material file Name
  //-----------------------
  fRun->SetMaterials("media_pnd.geo");
        
  // Create and add detectors
  //-------------------------
  FairModule *Cave= new PndCave("CAVE");
  Cave->SetGeometryFileName("pndcave.geo");
  fRun->AddModule(Cave); 
        
  FairModule *Magnet= new PndMagnet("MAGNET");
  Magnet->SetGeometryFileName("magnet.geo");
  fRun->AddModule(Magnet);
        
  FairModule *Pipe= new PndPipe("PIPE");
  Pipe->SetGeometryFileName("pipe.geo");
  fRun->AddModule(Pipe);

  PndEmc *Emc = new PndEmc("EMC",kTRUE);
  //new BwEndCap + FwEndCap
  //Emc->SetGeometryFileNameTriple("emc_module125.dat","emc_module3new.root","emc_module4_StraightGeo24.4.root"); //MapperVersion: 6
  Emc->SetGeometryVersion(15);
  Emc->SetStorageOfData(kFALSE);
  fRun->AddModule(Emc);
        
  // Create and Set Event Generator
  //-------------------------------
        
  FairPrimaryGenerator* primGen = new FairPrimaryGenerator();
  fRun->SetGenerator(primGen);
 
  TString part="gamma";
  cout <<" part = "<< part<<endl;
  TDatabasePDG *pdg = new TDatabasePDG();
  Int_t part_pid = pdg->GetParticle(part)->PdgCode();
  
  cout <<"pid == "<< part_pid<<endl;
  
  // Box Generator: 
  FairBoxGenerator* boxGen = new FairBoxGenerator(part_pid, 5); // 13 = muon; 1 = multipl. // 211 = pi+
  // first number: PDG particle code: 2nd number: particle multiplicity per event
        
  boxGen->SetPRange(1.,1.); // GeV/c
  // boxGen->SetPtRange(1.,1.); // GeV/c
  boxGen->SetPhiRange(0., 360.); // Azimuth angle range [degree]
  boxGen->SetThetaRange(0., 180.); // Polar angle in lab system range [degree]
  boxGen->SetXYZ(0., 0., 0.); // vertex coordinates [mm]
  primGen->AddGenerator(boxGen);  
        
  // fRun->SetStoreTraj(kTRUE); // to store particle trajectories 
  fRun->SetBeamMom(15); 
        //---------------------Create and Set the Field(s)---------- 
  PndMultiField *fField= new PndMultiField("FULL");
  fRun->SetField(fField);
        
  PndEmcHitProducer* emcHitProd = new PndEmcHitProducer();
  emcHitProd->SetStorageOfData(kFALSE);
  fRun->AddTask(emcHitProd);

  PndEmcHitsToWaveform* emcHitsToWaveform= new PndEmcHitsToWaveform(0);
  PndEmcWaveformToDigi* emcWaveformToDigi=new PndEmcWaveformToDigi(0);
  emcHitsToWaveform->SetStorageOfData(kFALSE);
  emcWaveformToDigi->SetStorageOfData(kFALSE);
  fRun->AddTask(emcHitsToWaveform);  // full digitization
  fRun->AddTask(emcWaveformToDigi);  // full digitization
 
  PndEmcMakeCluster* emcMakeCluster= new PndEmcMakeCluster(0);
  emcMakeCluster->SetStorageOfData(kFALSE);
  fRun->AddTask(emcMakeCluster);

  PndEmcHdrFiller* emcHdrFiller = new PndEmcHdrFiller();
  fRun->AddTask(emcHdrFiller); // ECM header

  PndEmcMakeBump* emcMakeBump= new PndEmcMakeBump();
  emcMakeBump->SetStorageOfData(kFALSE);
  fRun->AddTask(emcMakeBump);

  PndEmcMakeRecoHit* emcMakeRecoHit= new PndEmcMakeRecoHit();
  fRun->AddTask(emcMakeRecoHit);
   
   cout<< "emc_complete.C --> BEFORE corrections "<<endl;
  // ========================================================================
   PndEmcMakeCorr* emcCorrection= new PndEmcMakeCorr(2,transportModel);
   emcCorrection->SetStorageOfData(kFALSE);
   fRun->AddTask(emcCorrection);
  // ========================================================================


  fRun->Init();
  //PndEmcMapper *emcMap = PndEmcMapper::Instance(6);   
  PndMultiFieldPar* Par = (PndMultiFieldPar*) rtdb->getContainer("PndMultiFieldPar");
  if (fField) {  Par->SetParameters(fField); }
  Par->setInputVersion(fRun->GetRunId(),1);
  Par->setChanged();
        
  rtdb->saveOutput();
  //    rtdb->print();
                
  // Transport nEvents
  // -----------------
  cout <<" nEvents == "<<nEvents<<endl;
  fRun->Run(nEvents);
  cout <<" ***done == "<<endl;
  timer.Stop();
        
  Double_t rtime = timer.RealTime();
  Double_t ctime = timer.CpuTime();
  printf("RealTime=%f seconds, CpuTime=%f seconds\n",rtime,ctime);
        
}