//////////////////////////////////////////////////////////////////////////////////// // // _____ Source Position macro_____ // // Take as input root files with hillas parameters and recompute the ones depending // on the source position for a new (input) position, optionally rotating it in an // event by event basis. Output is a file with recomputed hillas parameters // // Ester Aliu // Oscar Blanch // Javier Rico //////////////////////////////////////////////////////////////////////////////////// #include #include #include "TString.h" #include "MHillasSrcCalc.h" #include "MHillasSrc.h" #include "MSrcRotate.h" #include "MSrcPosFromFile.h" #include "MSrcTranslate.h" #include "MParList.h" #include "MTaskList.h" #include "MHillas.h" #include "MReadTree.h" #include "MEvtLoop.h" #include "MLog.h" #include "MArgs.h" #include "MWriteRootFile.h" #include "MTime.h" #include "MControlPlots.h" #include "MIslands.h" using namespace std; Bool_t readDatacards(TString& filename); void srcPos(); //----------------------------------------------------------------------------- // declaration of variables read from datacards //----------------------------------------------------------------------------- TString inputFile; TString offFile; TString outputFile; TString offOutputFile; ULong_t nmaxevents=999999999; Float_t xsrcpos=0.; Float_t ysrcpos=0.; Double_t mjdpos=0; Bool_t kRotate=0; Bool_t kSrcPolicy=kFALSE; Double_t fRA= -1.; Double_t fDEC= -1.; TString srcFile; TString controlplotsfilename="controlplots.ps"; //----------------------------------------------------------------------------- // constants //----------------------------------------------------------------------------- const TString defaultcard="srcpos.datacard"; const Float_t conver = 189./0.6; // conversion factor degrees to mm //----------------------------------------------------------------------------- static void Usage() { gLog <" << endl << endl; gLog << " : datacards file name (dafault " << defaultcard <<")" << endl; gLog << " -?/-h: This help" << endl << endl; } //----------------------------------------------------------------------------- int main(int argc, char **argv) { // evaluate arguments MArgs arg(argc, argv); if (arg.HasOnly("-?") || arg.HasOnly("-h") || arg.HasOnly("--help")) { Usage(); return -1; } TString datacard = arg.GetArgumentStr(0); if(!datacard.Length()) datacard = defaultcard; if(!readDatacards(datacard)) { cout << "Error reading datacards. Stoping" << endl; return -1; } srcPos(); } //----------------------------------------------------------------------------- void srcPos() { // variable declaration Float_t xpos=xsrcpos*conver; // [mm] Float_t ypos=ysrcpos*conver; // [mm] // containers MParList plist; MTaskList tlist; MIslands islands; // include containers in parameter list plist.AddToList(&tlist); plist.AddToList(&islands); // tasks MReadTree read("Parameters", inputFile); read.DisableAutoScheme(); MSrcTranslate srctranslate; MSrcRotate srcrotate; MSrcPosFromFile srcfromfile; if(srcFile.Length()) { srcfromfile.SetInputFileName(srcFile); srcfromfile.SetMode(MSrcPlace::kOn); } else { srctranslate.SetTranslation(xpos,ypos); srctranslate.SetRelativeTranslation(kSrcPolicy); srctranslate.SetCreateHisto(kFALSE); srcrotate.SetRAandDECandRefMJD(fRA,fDEC,mjdpos); srcrotate.SetMode(MSrcPlace::kOn); } MHillasSrcCalc csrc1; MWriteRootFile write(outputFile,"RECREATE"); write.AddContainer("MHillas" , "Parameters"); write.AddContainer("MHillasSrc" , "Parameters"); write.AddContainer("MHillasExt" , "Parameters"); write.AddContainer("MNewImagePar" , "Parameters"); write.AddContainer("MRawEvtHeader" , "Parameters"); write.AddContainer("MRawRunHeader" , "Parameters"); write.AddContainer("MTime" , "Parameters"); write.AddContainer("MConcentration" , "Parameters"); write.AddContainer("MSrcPosCam" , "Parameters"); write.AddContainer("MIslands" , "Parameters"); //write.AddContainer("MIslands2" , "Parameters"); MControlPlots controlplots(controlplotsfilename); controlplots.SetProduceFile(kFALSE); // include tasks in task list tlist.AddToList(&read); if(srcFile.Length()) tlist.AddToList(&srcfromfile); else { tlist.AddToList(&srctranslate); if(kRotate) tlist.AddToList(&srcrotate); } tlist.AddToList(&csrc1); tlist.AddToList(&controlplots); tlist.AddToList(&write); // Eventloop MEvtLoop evtloop; evtloop.SetParList(&plist); if (!evtloop.Eventloop(nmaxevents)) return; tlist.PrintStatistics(); // do off-data if input file was specified if(!offFile.Length()) return; MReadTree read2("Parameters", offFile); read2.DisableAutoScheme(); tlist.AddToListBefore(&read2, &read, "All"); tlist.RemoveFromList(&read); MWriteRootFile write2(offOutputFile,"RECREATE"); write2.AddContainer("MHillas" , "Parameters"); write2.AddContainer("MHillasSrc" , "Parameters"); write2.AddContainer("MHillasExt" , "Parameters"); write2.AddContainer("MNewImagePar" , "Parameters"); write2.AddContainer("MRawEvtHeader" , "Parameters"); write2.AddContainer("MRawRunHeader" , "Parameters"); write2.AddContainer("MTime" , "Parameters"); write2.AddContainer("MConcentration" , "Parameters"); write2.AddContainer("MSrcPosCam" , "Parameters"); write2.AddContainer("MIslands" , "Parameters"); // write2.AddContainer("MIslands2" , "Parameters"); tlist.AddToListBefore(&write2,&write,"All"); tlist.RemoveFromList(&write); if(srcFile.Length()) srcfromfile.SetMode(MSrcPlace::kOff); else srcrotate.SetMode(MSrcPlace::kOff); controlplots.SetMode(MControlPlots::kOff); controlplots.SetProduceFile(kTRUE); if (!evtloop.Eventloop(nmaxevents)) return; tlist.PrintStatistics(); } //----------------------------------------------------------------------- Bool_t readDatacards(TString& filename) { ifstream ifun(filename.Data()); if(!ifun) { cout << "File " << filename << " not found" << endl; return kFALSE; } TString word; while(ifun >> word) { // skip comments if(word[0]=='/' && word[1]=='/') { while(ifun.get()!='\n'); // skip line continue; } // number of events if(strcmp(word.Data(),"NEVENTS")==0) ifun >> nmaxevents; // input file name if(strcmp(word.Data(),"INPUTFILES")==0) { if(inputFile.Length()) cout << "readDataCards Warning: overriding on-data file name" << endl; ifun >> inputFile; } // off-data file name if(strcmp(word.Data(),"OFFFILES")==0) { if(offFile.Length()) cout << "readDataCards Warning: overriding off-data file name" << endl; ifun >> offFile; } // output file name if(strcmp(word.Data(),"OUTFILE")==0) { if(outputFile.Length()) cout << "readDataCards Warning: overriding output file name" << endl; ifun >> outputFile; } // output file name (off data) if(strcmp(word.Data(),"OFFOUTFILE")==0) { if(offOutputFile.Length()) cout << "readDataCards Warning: overriding output file name for off data" << endl; ifun >> offOutputFile; } // source position input file if(strcmp(word.Data(),"SRCFILE")==0) { if(srcFile.Length()) cout << "readDataCards Warning: overriding source-position file name" << endl; ifun >> srcFile; } // source position if(strcmp(word.Data(),"SRCPOS")==0) { ifun >> xsrcpos; ifun >> ysrcpos; ifun >> mjdpos; } // source celestial coordinates if(strcmp(word.Data(),"SRCCOORDS")==0) { ifun >> fRA; ifun >> fDEC; } // source movement policy flag if(strcmp(word.Data(),"SRCABS")==0) ifun >> kSrcPolicy; // rotation flag if(strcmp(word.Data(),"ROTFLAG")==0) ifun >> kRotate; } // check compulsory values if(!inputFile.Length()) { cout << "No on-data file name specified" << endl; return kFALSE; } if(!outputFile.Length()) { cout << "No output file name specified" << endl; return kFALSE; } if(offFile.Length() && !offOutputFile.Length()) { cout << "No output file name specified for off data" << endl; return kFALSE; } if(xsrcpos==0 && ysrcpos==0) { cout << "Source position is center of the camera (as in input file)" << endl; return kFALSE; } MTime thetime(mjdpos); // Dump read values cout << "************************************************" << endl; cout << "* Datacards read from file " << filename << endl; cout << "************************************************" << endl; cout << "Maximum number of input events: " << nmaxevents << endl; cout << "Input file name(s): " << inputFile << endl; if(offFile.Length()) cout << "Input file name(s) for off data: " << offFile << endl; cout << "Output file name: " << outputFile << endl; if(offFile.Length()) cout << "Output file name for off data: " << offOutputFile << endl; if(srcFile.Length()) cout << "Reading source position from file " << srcFile << endl; else { cout << "Source position (degrees) X=" << xsrcpos << ", Y="<