runLumiPixel0SimBox.C

Go to the documentation of this file.

// Panda FullSim macro

//void runLumi0SimBox(const int nEvents=10, const double mom=15, TString storePath="tmpOutput", const int verboseLevel=0, const int particle=-211)
int runLumiPixel0SimBox(const int nEvents = 100, const int startEv = 0, TString storePath = "tmpOutput",
    const int verboseLevel = 0, const int particle = -2212, double mom = 15, const int trkNum = 1,
    const int seed = 0, const double dP = 0, TString geometryFile = "") {

        // ///PROOF lite
        // TProof::Open("");

        if (geometryFile == "") {
                geometryFile = "Luminosity-Detector.root";
        }

        gRandom->SetSeed(seed);
        //gRandom->SetSeed(0);
        TStopwatch timer;
        timer.Start();
        gDebug = 0;
        mom += dP;
        cout << "We start run for beam Mom = " << mom << endl;
        //output1
        TString simOutput = storePath + "/Lumi_MC_";
        simOutput += startEv;
        simOutput += ".root";
        TString parOutput = storePath + "/Lumi_Params_";
        parOutput += startEv;
        parOutput += ".root";

        // //Load basic libraries
        // gROOT->Macro("$VMCWORKDIR/gconfig/rootlogon.C");
        // gSystem->Load("libSds");
        // gSystem->Load("libLmd");
        FairRunSim *fRun = new FairRunSim();
        cout << "All libraries succsesfully loaded!" << endl;

        //set the MC version used
        fRun->SetName("TGeant4");
        //fRun->SetName("TGeant3");  //GEANE uses GEANT3!

        fRun->SetOutputFile(simOutput);

        //set material
        fRun->SetMaterials("media_pnd.geo");

        //create and add detectors
        FairModule *Cave = new PndCave("CAVE");
        //Cave->SetGeometryFileName("../macro/lmd/pndcaveVAC.geo");
        Cave->SetGeometryFileName("pndcave.geo");
        fRun->AddModule(Cave);

        FairModule *Pipe = new PndPipe("PIPE");
        //  Pipe->SetGeometryFileName("../macro/lmd/geo/beampipe_201303.root");
        Pipe->SetGeometryFileName("beampipe_201309.root");
        fRun->AddModule(Pipe);

        FairModule *Magnet = new PndMagnet("MAGNET");
        Magnet->SetGeometryFileName("FullSolenoid_V842.root");
        fRun->AddModule(Magnet);

        FairModule *Dipole = new PndMagnet("MAGNET");
        Dipole->SetGeometryFileName("dipole.geo");
        fRun->AddModule(Dipole);

        PndLmdDetector *Lum = new PndLmdDetector("LUM", kTRUE);
        Lum->SetExclusiveSensorType("LumActive");  //ignore MVD
        //  Lum->SetGeometryFileName("../macro/lmd/geo/Test-Dipol-Design.root"); //sensors with trap shape
        //  Lum->SetGeometryFileName("../macro/lmd/geo/HV_MAPS-Design-29052013.root"); // LMD including box etc
        Lum->SetGeometryFileName(geometryFile);
        //Lum->SetGeometryFileName("../macro/lmd/geo/HV_MAPS-Design-SensorsOnly.root"); // LMD, sensors only
        Lum->SetVerboseLevel(verboseLevel);

        fRun->AddModule(Lum);

        //particle generator
        FairPrimaryGenerator* primGen = new FairPrimaryGenerator();

        fRun->SetGenerator(primGen);

        // Box Generator
        FairBoxGenerator *fBox = new FairBoxGenerator(particle, trkNum);
        fBox->SetPRange(mom, mom);
        //  fBox->SetThetaRange(0.52,0.63); // 9 ... 11 mrad
        //  fBox->SetThetaRange(0.12,0.7); // 2... 12 mrad
        //  fBox->SetThetaRange(1.3,1.4); // mrad outside of the detector geometry by purpose
        fBox->SetThetaRange(0.13, 0.65);  // 2... 11 mrad
        //fBox->SetThetaRange(0.12,0.65); // 2... 11 mrad
        //fBox->SetThetaRange(0.229183, 0.458366); //4 ... 8 mrad
        //fBox->SetThetaRange(0.229183,0.31512);//4..5.5 mrad
        //fBox->SetThetaRange(0.229,0.229);//4..mrad
        //fBox->SetThetaRange(0.458366, 0.458366); //!!! 8 mrad
        //fBox->SetThetaRange(0., 45.);//TEST
        //fBox->SetPhiRange(90,90.);
        fBox->SetPhiRange(0, 360.);
        //fBox->SetPhiRange(0.5,359.5); //FOR missed track check
        //fBox->SetPhiRange(0.,20.);//!!! TEST
        //fBox->SetThetaRange(0.344,0.344); //!!! TEST ~ 6 mrad
        //fBox->SetPhiRange(45,45);//TEST
        primGen->AddGenerator(fBox);

        //reading the new field map in the old format
        fRun->SetBeamMom(mom);

        PndMultiField *fField = new PndMultiField("AUTO");

        fRun->SetField(fField);

        if (nEvents < 101) fRun->SetStoreTraj(kTRUE);  // toggle this for use with EVE
        else fRun->SetStoreTraj(kFALSE);

        bool misalignedGeometry = true;
        if (misalignedGeometry) {
                string misMatricesFilePath = "geo/misalignMatrices-SensorsOnly-100.root";
                // check if file exists, if true, try to read it
                TFile *misalignmentMatrixRootfile = new TFile(misMatricesFilePath.c_str(), "READ");
                if (misalignmentMatrixRootfile->IsOpen()) {
                        printf("File opened successfully\n");
                        std::map < std::string, TGeoHMatrix > *matrices;
                        gDirectory->GetObject("PndLmdMisalignMatrices", matrices);
                        misalignmentMatrixRootfile->Close();
                        cout << matrices->size() << " matrices successfully read from file.";
                        Lum->SetMisalignmentMatrices(*matrices);
                        cout << "matrix set!\n";
                }

                // if not, fail violently
                else {
                        cerr << "WARNING. I was instructed to use misaligned geometry,\n";
                        cerr << "but no misaligned matrices could be found in " << misMatricesFilePath << "\n";
                        return 1;
                }
        }

        fRun->Init();

        // Fill the Parameter containers for this run
        //-------------------------------------------
        FairRuntimeDb *rtdb = fRun->GetRuntimeDb();
        Bool_t kParameterMerged = kTRUE;
        FairParRootFileIo* output = new FairParRootFileIo(kParameterMerged);
        output->open(parOutput.Data(), "RECREATE");
        rtdb->setOutput(output);
        PndMultiFieldPar* Par = (PndMultiFieldPar*) rtdb->getContainer("PndMultiFieldPar");
        if (fField) {
                Par->SetParameters(fField);
        }
        Par->setInputVersion(fRun->GetRunId(), 1);
        Par->setChanged();

        // Transport nEvents
        // -----------------

        fRun->Run(nEvents);

        rtdb->saveOutput();
        rtdb->print();

        timer.Stop();
        Double_t rtime = timer.RealTime();
        Double_t ctime = timer.CpuTime();
        printf("RealTime=%f seconds, CpuTime=%f seconds\n", rtime, ctime);

        FairLogger *logger = FairLogger::GetLogger();

        // log to screen and to file
        logger->SetLogToScreen(kTRUE);
        logger->SetLogToFile(kFALSE);

        logger->SetLogVerbosityLevel("LOW");

        // Set different levels of verbosity. In the example everything >=INFO goes to the
        // file and everything >= ERROR is printed on the screen
        // LogLevels are (FATAL, ERROR, WARNING, INFO, DEBUG, DEBUG1, DEBUG2, DEBUG3, DEBUG4)
        logger->SetLogScreenLevel("INFO");  //Only FATAL and ERROR to screen

        // temporary fix to avoid double frees at the destruction of te program for pandaroot/fairroot with root6
        gGeoManager->GetListOfVolumes()->Delete();
        gGeoManager->GetListOfShapes()->Delete();
        delete gGeoManager;

        return 0;
}