estarlight/src/e_starlightStandalone.cpp

0001 //////////////////////////////////////////////////////////////////////////
0002 //
0003 //    Copyright 2010
0004 //
0005 //    This file is part of starlight.
0006 //
0007 //    starlight is free software: you can redistribute it and/or modify
0008 //    it under the terms of the GNU General Public License as published by
0009 //    the Free Software Foundation, either version 3 of the License, or
0010 //    (at your option) any later version.
0011 //
0012 //    starlight is distributed in the hope that it will be useful,
0013 //    but WITHOUT ANY WARRANTY; without even the implied warranty of
0014 //    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
0015 //    GNU General Public License for more details.
0016 //
0017 //    You should have received a copy of the GNU General Public License
0018 //    along with starlight. If not, see <http://www.gnu.org/licenses/>.
0019 //
0020 ///////////////////////////////////////////////////////////////////////////
0021 //
0022 // File and Version Information:
0023 // $Rev:: 270                         $: revision of last commit
0024 // $Author:: jnystrand                $: author of last commit
0025 // $Date:: 2016-07-08 16:31:51 +0100 #$: date of last commit
0026 //
0027 // Description:
0028 //
0029 //
0030 //
0031 ///////////////////////////////////////////////////////////////////////////
0032
0033
0034 #include <iostream>
0035 #include <chrono>
0036 #include "reportingUtils.h"
0037 #include "e_starlight.h"
0038 #include "inputParameters.h"
0039 #include "eventfilewriter.h"
0040 #include "e_starlightStandalone.h"
0041
0042
0043 using namespace std;
0044
0045
0046 e_starlightStandalone::e_starlightStandalone()
0047     :   _configFileName   ("slight.in"),
0048         _starlight        (0),
0049         _nmbEventsTot     (1),
0050         _nmbEventsPerFile (_nmbEventsTot)
0051 { }
0052
0053
0054 e_starlightStandalone::~e_starlightStandalone()
0055 { }
0056
0057
0058 bool
0059 e_starlightStandalone::init()
0060 {
0061     _inputParameters = new inputParameters();
0062     // read input parameters from config file
0063     _inputParameters->configureFromFile(_configFileName);
0064     if (!_inputParameters->init()) {
0065         printWarn << "problems initializing input parameters. cannot initialize starlight." << endl;
0066         return false;
0067     }
0068
0069     // copy input file to one with baseFileName naming scheme
0070         std::string inputCopyName, _baseFileName;
0071         _baseFileName = _inputParameters->baseFileName();
0072     if (_baseFileName != "slight") {
0073          inputCopyName = _baseFileName +".in";
0074
0075         ofstream inputCopyFile;
0076     inputCopyFile.open(inputCopyName.c_str());
0077
0078         std::ifstream infile(_configFileName.c_str());
0079         if ((!infile) || (!infile.good()))
0080         {
0081           return -1;
0082         }
0083
0084         int lineSize = 256;
0085         char tmp[lineSize];
0086         while (!infile.getline(tmp, lineSize).eof())
0087          {
0088          cout << tmp << endl;
0089          inputCopyFile << tmp << endl;
0090          }
0091         inputCopyFile.close();
0092     }
0093
0094     // get the number of events
0095     // for now we write everything to one file
0096     _nmbEventsTot     = _inputParameters->nmbEvents();
0097     _nmbEventsPerFile = _nmbEventsTot;
0098
0099     // create the starlight object
0100     _starlight = new e_starlight();
0101
0102         // give starlight the input parameters
0103         _starlight->setInputParameters(_inputParameters);
0104
0105     // initialize starlight
0106     return _starlight->init();
0107 }
0108
0109
0110 bool
0111 e_starlightStandalone::run()
0112 {
0113     if (!_starlight) {
0114         printWarn << "null pointer to starlight object. make sure that init() was called. "
0115                   << "cannot generate events." << endl;
0116         return false;
0117     }
0118
0119     // open output file
0120     eventFileWriter fileWriter;
0121     eventFileWriter fileWriterLUND;
0122     fileWriter.writeFullPythiaInfo(_inputParameters->pythiaFullEventRecord());
0123     fileWriter.writeFullHepMC3Info(_inputParameters->hepmc3FullEventRecord());
0124         _baseFileName = _inputParameters->baseFileName();
0125         _eventDataFileName = _baseFileName +".out";
0126     fileWriter.open(_eventDataFileName);
0127     //
0128     fileWriter.writeInit(*_inputParameters);
0129     //
0130     if(_inputParameters->lundFullEventRecord()){
0131         fileWriterLUND.open("slight_LUND.txt");
0132         //
0133         fileWriterLUND.writeInitLUND(*_inputParameters);
0134         //
0135     }
0136     printInfo << "generating events:" << endl;
0137     unsigned int nmbEvents = 0;
0138     std::chrono::steady_clock::time_point begin= std::chrono::steady_clock::now();
0139     while (nmbEvents < _nmbEventsTot) {
0140         for (unsigned int iEvent = 0; (iEvent < _nmbEventsPerFile) && (nmbEvents < _nmbEventsTot);
0141              ++iEvent, ++nmbEvents) {
0142             progressIndicator(iEvent, _nmbEventsTot, true, 4);
0143             eXEvent event = _starlight->produceEvent();
0144             // Boost event from back to lab reference frame
0145             boostEvent(event);
0146             reflectEvent(event);
0147             fileWriter.writeEvent(event, iEvent);
0148             if(_inputParameters->lundFullEventRecord()){
0149                 fileWriterLUND.writeEventLUND(event, iEvent);
0150             }
0151         }
0152     }
0153     std::chrono::steady_clock::time_point end= std::chrono::steady_clock::now();
0154     float running_total = 1E-3*std::chrono::duration_cast<std::chrono::milliseconds>(end - begin).count();
0155     cout<<"Total time "<<running_total<<" s ("<<1E3*running_total/nmbEvents<<" ms/ev)"<<endl;
0156     fileWriter.close();
0157     if(_inputParameters->lundFullEventRecord()){
0158         fileWriterLUND.close();
0159     }
0160
0161     if( _starlight->nmbAttempts() == 0 )return true;
0162
0163     double _branchingRatio = _inputParameters->inputBranchingRatio();
0164     printInfo << "number of attempts = " << _starlight->nmbAttempts() << ", "
0165               << "number of accepted events = " << _starlight->nmbAccepted() << endl;
0166         double selectedCrossSection =
0167       ((double)_starlight->nmbAccepted()/_starlight->nmbAttempts())*_branchingRatio*_starlight->getTotalCrossSection();
0168     if (selectedCrossSection > 1.){
0169       cout<< " The cross section of the generated sample is "<<selectedCrossSection<<" barn."<<endl;
0170     } else if (1.E3*selectedCrossSection > 1.){
0171       cout<< " The cross section of the generated sample is "<<1.E3*selectedCrossSection<<" mb."<<endl;
0172         } else if (1.E6*selectedCrossSection > 1.){
0173       cout<< " The cross section of the generated sample is "<<1.E6*selectedCrossSection<<" microbarn."<<endl;
0174         } else if (1.E9*selectedCrossSection > 1.){
0175       cout<< " The cross section of the generated sample is "<<1.E9*selectedCrossSection<<" nanobarn."<<endl;
0176         } else if (1.E12*selectedCrossSection > 1.){
0177       cout<< " The cross section of the generated sample is "<<1.E12*selectedCrossSection<<" picobarn."<<endl;
0178         } else {
0179       cout<< " The cross section of the generated sample is " <<1.E15*selectedCrossSection<<" femtob."<<endl;
0180         }
0181
0182     return true;
0183 }
0184 void e_starlightStandalone::boostEvent(eXEvent &event)
0185 {
0186
0187    // Calculate CMS boost
0188    double rap2 = -acosh(_inputParameters->targetBeamLorentzGamma());
0189    double boost = _inputParameters->rap_CM();
0190    event.boost(boost, rap2); //  Boost back to laboratory reference frame. Electron initially in target frame and V.M. in CMS
0191                              //  Assuming electron is electronBeam
0192
0193
0194 }
0195
0196
0197 void e_starlightStandalone::reflectEvent(eXEvent &event)
0198 {
0199    event.reflect(); //Change all z to -z
0200 }