source: trunk/MagicSoft/Mars/mtemp/mifae/programs/srcPos.cc@ 4074

Last change on this file since 4074 was 4074, checked in by rico, 21 years ago
*** empty log message ***
File size: 9.1 KB
Line 
1////////////////////////////////////////////////////////////////////////////////////
2//
3// _____ Source Position macro_____
4//
5// Take as input root files with hillas parameters and recompute the ones depending
6// on the source position for a new (input) position, optionally rotating it in an
7// event by event basis. Output is a file with recomputed hillas parameters
8//
9// Ester Aliu <aliu@ifae.es>
10// Oscar Blanch <blanch@ifae.es>
11// Javier Rico <jrico@ifae.es>
12////////////////////////////////////////////////////////////////////////////////////
13
14#include <fstream>
15#include <iostream>
16
17#include "TString.h"
18
19#include "MHillasSrcCalc.h"
20#include "MHillasSrc.h"
21#include "MSrcRotate.h"
22#include "MSrcPosFromFile.h"
23#include "MSrcTranslate.h"
24#include "MParList.h"
25#include "MTaskList.h"
26#include "MHillas.h"
27#include "MReadTree.h"
28#include "MEvtLoop.h"
29#include "MLog.h"
30#include "MArgs.h"
31#include "MWriteRootFile.h"
32#include "MTime.h"
33
34using namespace std;
35
36Bool_t readDatacards(TString& filename);
37void srcPos();
38
39//-----------------------------------------------------------------------------
40// declaration of variables read from datacards
41//-----------------------------------------------------------------------------
42
43TString inputFile;
44TString offFile;
45TString outputFile;
46TString offOutputFile;
47ULong_t nmaxevents=999999999;
48Float_t xsrcpos=0.;
49Float_t ysrcpos=0.;
50Double_t mjdpos=0;
51Bool_t kRotate=0;
52Bool_t kSrcPolicy=kFALSE;
53Double_t fRA= -1.;
54Double_t fDEC= -1.;
55TString srcFile;
56
57//-----------------------------------------------------------------------------
58// constants
59//-----------------------------------------------------------------------------
60
61const TString defaultcard="srcpos.datacard";
62const Float_t conver = 189./0.6; // conversion factor degrees to mm
63
64//-----------------------------------------------------------------------------
65
66static void Usage()
67{
68 gLog <<endl;
69 gLog << "Usage is:" << endl;
70 gLog << " srcPos [-h] [-?] <datacards>" << endl << endl;
71 gLog << " <datacards>: datacards file name (dafault " << defaultcard <<")" << endl;
72 gLog << " -?/-h: This help" << endl << endl;
73}
74
75//-----------------------------------------------------------------------------
76int main(int argc, char **argv)
77{
78 // evaluate arguments
79 MArgs arg(argc, argv);
80 if (arg.HasOnly("-?") || arg.HasOnly("-h") || arg.HasOnly("--help"))
81 {
82 Usage();
83 return -1;
84 }
85
86 TString datacard = arg.GetArgumentStr(0);
87 if(!datacard.Length())
88 datacard = defaultcard;
89
90 if(!readDatacards(datacard))
91 {
92 cout << "Error reading datacards. Stoping" << endl;
93 return -1;
94 }
95 srcPos();
96}
97
98//-----------------------------------------------------------------------------
99void srcPos()
100{
101 // variable declaration
102 Float_t xpos=xsrcpos*conver; // [mm]
103 Float_t ypos=ysrcpos*conver; // [mm]
104
105 // containers
106 MParList plist;
107 MTaskList tlist;
108
109 // include containers in parameter list
110 plist.AddToList(&tlist);
111
112 // tasks
113 MReadTree read("Parameters", inputFile);
114 read.DisableAutoScheme();
115
116 MSrcTranslate srctranslate;
117 MSrcRotate srcrotate;
118 MSrcPosFromFile srcfromfile;
119
120 if(srcFile.Length())
121 {
122 srcfromfile.SetInputFileName(srcFile);
123 srcfromfile.SetMode(MSrcPlace::kOn);
124 }
125 else
126 {
127 srctranslate.SetTranslation(xpos,ypos);
128 srctranslate.SetRelativeTranslation(kSrcPolicy);
129 srcrotate.SetRAandDECandRefMJD(fRA,fDEC,mjdpos);
130 srcrotate.SetMode(MSrcPlace::kOn);
131 }
132
133 MHillasSrcCalc csrc1;
134
135 MWriteRootFile write(outputFile,"RECREATE");
136 write.AddContainer("MHillas" , "Parameters");
137 write.AddContainer("MHillasSrc" , "Parameters");
138 write.AddContainer("MHillasExt" , "Parameters");
139 write.AddContainer("MNewImagePar" , "Parameters");
140 write.AddContainer("MRawEvtHeader" , "Parameters");
141 write.AddContainer("MRawRunHeader" , "Parameters");
142 write.AddContainer("MConcentration" , "Parameters");
143 write.AddContainer("MSrcPosCam" , "Parameters");
144
145 // include tasks in task list
146 tlist.AddToList(&read);
147 if(srcFile.Length())
148 tlist.AddToList(&srcfromfile);
149 else
150 {
151 tlist.AddToList(&srctranslate);
152 if(kRotate)
153 tlist.AddToList(&srcrotate);
154 }
155 tlist.AddToList(&csrc1);
156 tlist.AddToList(&write);
157
158 // Eventloop
159 MEvtLoop evtloop;
160 evtloop.SetParList(&plist);
161 if (!evtloop.Eventloop(nmaxevents))
162 return;
163
164 tlist.PrintStatistics();
165
166
167 // do off-data if input file was specified
168 if(!offFile.Length())
169 return;
170
171 MReadTree read2("Parameters", offFile);
172 read2.DisableAutoScheme();
173 tlist.AddToListBefore(&read2, &read, "All");
174 tlist.RemoveFromList(&read);
175
176 MWriteRootFile write2(offOutputFile,"RECREATE");
177 write2.AddContainer("MHillas" , "Parameters");
178 write2.AddContainer("MHillasSrc" , "Parameters");
179 write2.AddContainer("MHillasExt" , "Parameters");
180 write2.AddContainer("MNewImagePar" , "Parameters");
181 write2.AddContainer("MRawEvtHeader" , "Parameters");
182 write2.AddContainer("MRawRunHeader" , "Parameters");
183 write2.AddContainer("MConcentration" , "Parameters");
184 write2.AddContainer("MSrcPosCam" , "Parameters");
185 tlist.AddToListBefore(&write2,&write,"All");
186 tlist.RemoveFromList(&write);
187
188 if(srcFile.Length())
189 srcfromfile.SetMode(MSrcPlace::kOff);
190 else
191 srcrotate.SetMode(MSrcPlace::kOff);
192
193 if (!evtloop.Eventloop(nmaxevents))
194 return;
195
196 tlist.PrintStatistics();
197
198}
199//-----------------------------------------------------------------------
200
201Bool_t readDatacards(TString& filename)
202{
203 ifstream ifun(filename.Data());
204 if(!ifun)
205 {
206 cout << "File " << filename << " not found" << endl;
207 return kFALSE;
208 }
209
210 TString word;
211
212 while(ifun >> word)
213 {
214 // skip comments
215 if(word[0]=='/' && word[1]=='/')
216 {
217 while(ifun.get()!='\n'); // skip line
218 continue;
219 }
220
221 // number of events
222 if(strcmp(word.Data(),"NEVENTS")==0)
223 ifun >> nmaxevents;
224
225 // input file name
226 if(strcmp(word.Data(),"INPUTFILES")==0)
227 {
228 if(inputFile.Length())
229 cout << "readDataCards Warning: overriding on-data file name" << endl;
230 ifun >> inputFile;
231 }
232
233 // off-data file name
234 if(strcmp(word.Data(),"OFFFILES")==0)
235 {
236 if(offFile.Length())
237 cout << "readDataCards Warning: overriding off-data file name" << endl;
238 ifun >> offFile;
239 }
240
241 // output file name
242 if(strcmp(word.Data(),"OUTFILE")==0)
243 {
244 if(outputFile.Length())
245 cout << "readDataCards Warning: overriding output file name" << endl;
246 ifun >> outputFile;
247 }
248
249 // output file name (off data)
250 if(strcmp(word.Data(),"OFFOUTFILE")==0)
251 {
252 if(offOutputFile.Length())
253 cout << "readDataCards Warning: overriding output file name for off data" << endl;
254 ifun >> offOutputFile;
255 }
256
257 // source position input file
258 if(strcmp(word.Data(),"SRCFILE")==0)
259 {
260 if(srcFile.Length())
261 cout << "readDataCards Warning: overriding source-position file name" << endl;
262 ifun >> srcFile;
263 }
264
265 // source position
266 if(strcmp(word.Data(),"SRCPOS")==0)
267 {
268 ifun >> xsrcpos;
269 ifun >> ysrcpos;
270 ifun >> mjdpos;
271 }
272
273 // source celestial coordinates
274 if(strcmp(word.Data(),"SRCCOORDS")==0)
275 {
276 ifun >> fRA;
277 ifun >> fDEC;
278 }
279
280 // source movement policy flag
281 if(strcmp(word.Data(),"SRCABS")==0)
282 ifun >> kSrcPolicy;
283
284 // rotation flag
285 if(strcmp(word.Data(),"ROTFLAG")==0)
286 ifun >> kRotate;
287 }
288
289 // check compulsory values
290 if(!inputFile.Length())
291 {
292 cout << "No on-data file name specified" << endl;
293 return kFALSE;
294 }
295
296 if(!outputFile.Length())
297 {
298 cout << "No output file name specified" << endl;
299 return kFALSE;
300 }
301
302 if(offFile.Length() && !offOutputFile.Length())
303 {
304 cout << "No output file name specified for off data" << endl;
305 return kFALSE;
306 }
307
308 if(xsrcpos==0 && ysrcpos==0)
309 {
310 cout << "Source position is center of the camera (as in input file)" << endl;
311 return kFALSE;
312 }
313
314 MTime thetime(mjdpos);
315
316 // Dump read values
317 cout << "************************************************" << endl;
318 cout << "* Datacards read from file " << filename << endl;
319 cout << "************************************************" << endl;
320 cout << "Maximum number of input events: " << nmaxevents << endl;
321 cout << "Input file name(s): " << inputFile << endl;
322 if(offFile.Length())
323 cout << "Input file name(s) for off data: " << offFile << endl;
324 cout << "Output file name: " << outputFile << endl;
325 if(offFile.Length())
326 cout << "Output file name for off data: " << offOutputFile << endl;
327 if(srcFile.Length())
328 cout << "Reading source position from file " << srcFile << endl;
329 else
330 {
331 cout << "Source position (degrees) X=" << xsrcpos << ", Y="<<ysrcpos;
332 if(kSrcPolicy)
333 cout << " (RELATIVE TO INITIAL SOURCE POSITION)";
334 else
335 cout << " (ABSOLUTE POSITION IN THE CAMERA)";
336 cout << ", at " << thetime.GetSqlDateTime() << endl;
337 cout << "De-rotation flag " << kRotate << endl;
338 cout << "Source celestial coordiantes (rad): RA = " << fRA << ", DEC = " << fDEC << endl;
339 }
340 cout << "***********" << endl << endl;
341
342 return kTRUE;
343}
Note: See TracBrowser for help on using the repository browser.