source: trunk/MagicSoft/Mars/manalysis/MSourcePosfromStarPos.cc@ 3491

Last change on this file since 3491 was 3488, checked in by wittek, 22 years ago
*** empty log message ***
File size: 23.3 KB
Line 
1/* ======================================================================== *\
2!
3! *
4! * This file is part of MARS, the MAGIC Analysis and Reconstruction
5! * Software. It is distributed to you in the hope that it can be a useful
6! * and timesaving tool in analysing Data of imaging Cerenkov telescopes.
7! * It is distributed WITHOUT ANY WARRANTY.
8! *
9! * Permission to use, copy, modify and distribute this software and its
10! * documentation for any purpose is hereby granted without fee,
11! * provided that the above copyright notice appear in all copies and
12! * that both that copyright notice and this permission notice appear
13! * in supporting documentation. It is provided "as is" without express
14! * or implied warranty.
15! *
16!
17! Author(s): Wolfgang Wittek 02/2004 <mailto:wittek@mppmu.mpg.de>
18!
19! Copyright: MAGIC Software Development, 2000-2004
20!
21!
22\* ======================================================================== */
23
24/////////////////////////////////////////////////////////////////////////////
25//
26// MSourcePosfromStarPos
27//
28// This is a task which
29// - calculates the position of the source in the camera
30// from the position of a known star in the camera
31// - and puts the source position into the container MSrcPosCam
32//
33// Input :
34// ASCII file containing for each run
35// - the run number
36// - the direction (theta, phi) the telescope is pointing to in [deg]
37// - the position (xStar, yStar) of a known star in the camera in [mm]
38// - the error (dxStar, dyStar) of this position in [mm]
39//
40// Output Containers :
41// MSrcPosCam
42//
43/////////////////////////////////////////////////////////////////////////////
44#include <TList.h>
45#include <TSystem.h>
46#include <TMatrixD.h>
47
48#include <fstream>
49
50#include "MSourcePosfromStarPos.h"
51
52#include "MParList.h"
53#include "MRawRunHeader.h"
54#include "MGeomCam.h"
55#include "MSrcPosCam.h"
56#include "MMcEvt.hxx"
57
58#include "MLog.h"
59#include "MLogManip.h"
60#include "MObservatory.h"
61
62ClassImp(MSourcePosfromStarPos);
63
64using namespace std;
65
66// --------------------------------------------------------------------------
67//
68// Default constructor.
69//
70MSourcePosfromStarPos::MSourcePosfromStarPos(
71 const char *name, const char *title)
72 : fIn(NULL)
73{
74 fName = name ? name : "MSourcePosfromStarPos";
75 fTitle = title ? title : "Calculate source position from star position";
76
77 fFileNames = new TList;
78 fFileNames->SetOwner();
79
80 fRuns = 0;
81 fSize = 100;
82
83 fRunNr.Set(fSize);
84
85 fThetaTel.Set(fSize);
86 fPhiTel.Set(fSize);
87 fdThetaTel.Set(fSize);
88 fdPhiTel.Set(fSize);
89
90 fDecStar.Set(1);
91 fRaStar.Set(1);
92 fxStar.ResizeTo(1,fSize);
93 fyStar.ResizeTo(1,fSize);
94 fdxStar.ResizeTo(1,fSize);
95 fdyStar.ResizeTo(1,fSize);
96
97 fStars = 0;
98 fDecSource = 0.0;
99 fRaSource = 0.0;
100}
101
102// --------------------------------------------------------------------------
103//
104// Delete the filename list and the input stream if one exists.
105//
106MSourcePosfromStarPos::~MSourcePosfromStarPos()
107{
108 delete fFileNames;
109 if (fIn)
110 delete fIn;
111}
112
113
114// --------------------------------------------------------------------------
115//
116// Set the sky coordinates of the source and of the star
117//
118// Input :
119// declination in units of (Deg, Min, Sec)
120// right ascension in units of (Hour, Min, Sec)
121//
122
123void MSourcePosfromStarPos::SetSourceAndStarPosition(
124 TString nameSource,
125 Double_t decSourceDeg, Double_t decSourceMin, Double_t decSourceSec,
126 Double_t raSourceHour, Double_t raSourceMin, Double_t raSourceSec,
127 TString nameStar,
128 Double_t decStarDeg, Double_t decStarMin, Double_t decStarSec,
129 Double_t raStarHour, Double_t raStarMin, Double_t raStarSec )
130{
131 *fLog << "MSourcePosfromStarPos::SetSourceAndStarPosition :" << endl;
132 *fLog << " Source : " << nameSource << " " << decSourceDeg << ":"
133 << decSourceMin << ":" << decSourceSec << endl;
134 *fLog << " Star : " << nameStar << " " << decStarDeg << ":"
135 << decStarMin << ":" << decStarSec << endl;
136
137 // convert into radians
138 fDecSource = (decSourceDeg + decSourceMin/60.0 + decSourceSec/3600.0)
139 / kRad2Deg;
140 fRaSource = (raSourceHour + raSourceMin/60.0 + raSourceSec/3600.0)
141 * 360.0 / (24.0 * kRad2Deg);
142
143 fStars += 1;
144 fDecStar.Set(fStars);
145 fDecStar[fStars-1] = (decStarDeg + decStarMin/60.0 + decStarSec/3600.0)
146 / kRad2Deg;
147 fRaStar[fStars-1] = (raStarHour + raStarMin/60.0 + raStarSec/3600.0)
148 * 360.0 / (24.0 * kRad2Deg);
149
150 *fLog << all << "MSourcePosfromStarPos::SetSourceAndStarPosition; fDecSource, fRaSource, fDecStar, fRaStar were set to : [radians] "
151 << fDecSource << ", " << fRaSource << ", "
152 << fDecStar[fStars-1] << ", " << fRaStar[fStars-1] << endl;
153}
154
155// --------------------------------------------------------------------------
156//
157// Set the sky coordinates of another star
158//
159// Input :
160// declination in units of (Deg, Min, Sec)
161// right ascension in units of (Hour, Min, Sec)
162//
163
164void MSourcePosfromStarPos::AddStar(
165 TString nameStar,
166 Double_t decStarDeg, Double_t decStarMin, Double_t decStarSec,
167 Double_t raStarHour, Double_t raStarMin, Double_t raStarSec )
168{
169 *fLog << "MSourcePosfromStarPos::AddStar :" << endl;
170 *fLog << " Star : " << nameStar << " " << decStarDeg << ":"
171 << decStarMin << ":" << decStarSec << endl;
172
173 // convert into radians
174 fStars = fDecStar.GetSize() + 1;
175 fDecStar.Set(fStars);
176 fDecStar[fStars-1] = (decStarDeg + decStarMin/60.0 + decStarSec/3600.0)
177 / kRad2Deg;
178 fRaStar[fStars-1] = (raStarHour + raStarMin/60.0 + raStarSec/3600.0)
179 * 360.0 / (24.0 * kRad2Deg);
180
181 *fLog << all << "MSourcePosfromStarPos::AddStar; fDecStar, fRaStar were set to : [radians] "
182 << fDecStar[fStars-1] << ", " << fRaStar[fStars-1] << endl;
183}
184
185// --------------------------------------------------------------------------
186//
187//
188Int_t MSourcePosfromStarPos::PreProcess(MParList *pList)
189{
190 MGeomCam *geom = (MGeomCam*)pList->FindObject(AddSerialNumber("MGeomCam"));
191 if (!geom)
192 {
193 *fLog << err << "MSourcePosfromStarPos : MGeomCam (Camera Geometry) missing in Parameter List... aborting." << endl;
194 return kFALSE;
195 }
196 fMm2Deg = geom->GetConvMm2Deg();
197 // fDistCameraReflector is the distance of the camera from the reflector center in [mm]
198 fDistCameraReflector = kRad2Deg / fMm2Deg;
199 *fLog << all << "MSourcePosfromStarPos::PreProcess; fMm2Deg, fDistCameraReflector = "
200 << fMm2Deg << ", " << fDistCameraReflector << endl;
201
202 fObservatory = (MObservatory*)pList->FindObject(AddSerialNumber("MObservatory"));
203 if (!fObservatory)
204 {
205 *fLog << err << "MObservatory not found... aborting" << endl;
206 return kFALSE;
207 }
208
209 fRun = (MRawRunHeader*)pList->FindObject("MRawRunHeader");
210 if (!fRun)
211 {
212 *fLog << err << "MSourcePosfromStarPos::PreProcess; MRawRunHeader not found... aborting." << endl;
213 return kFALSE;
214 }
215
216
217 fMcEvt = (MMcEvt*)pList->FindCreateObj("MMcEvt");
218 if (!fMcEvt)
219 {
220 *fLog << err << "MSourcePosfromStarPos::PreProcess; MMcEvt not found... aborting." << endl;
221 return kFALSE;
222 }
223
224
225 fSrcPos = (MSrcPosCam*)pList->FindObject(AddSerialNumber("MSrcPosCam"));
226 if (!fSrcPos)
227 {
228 *fLog << err << "MSourcePosfromStarPos::PreProcess; MSrcPosCam not found... aborting" << endl;
229 return kFALSE;
230 }
231
232 //---------------------------------------------------------------------
233 // read all files and call ReadData() to read and store the information
234 //
235
236 *fLog << all << "---------------------------------" << endl;
237 while(1)
238 {
239 if (!OpenNextFile())
240 {
241 *fLog << "there is no more file to open" << endl;
242 break;
243 }
244
245 *fLog << "read data" << endl;
246 while (1)
247 {
248 if (fIn->eof())
249 {
250 *fLog << "eof encountered; open next file" << endl;
251
252 if (!OpenNextFile()) break;
253 }
254
255 // FIXME! Set InputStreamID
256
257 ReadData();
258 }
259 }
260
261 *fLog << "all data were read" << endl;
262 FixSize();
263
264 if (fDecSource == 0.0 || fRaSource == 0.0 || fStars == 0)
265 {
266 *fLog << err << "MSourcePosfromStarPos::PreProcess; there are no sky coordinates defined for the source or from stars; fStars = " << fStars
267 << endl;
268 }
269 *fLog << all << "---------------------------------" << endl;
270
271 //-------------------------------------------------------------
272
273 return kTRUE;
274}
275
276//=========================================================================
277//
278// SourcefromStar
279//
280// this routine calculates the position of a source (for example Crab) in the camera
281// from the position of a star (for example ZetaTauri) in the camera. The latter
282// position may have been determined by analysing the DC currents in the different
283// pixels.
284//
285// Input : thetaTel, phiTel the direction the telescope is pointing to,
286// in local coordinates
287// f the distance between camera and reflector
288// decStar, raStar the position of the star in sky coordinates
289// decSource, raSource the position of the source in sky coordinates
290// xStar, yStar the position of the star in the camera
291// dxStar, dyStar error of the position of the star in the camera
292//
293// Output : xSource, ySource the calculated position of the source in the camera
294// dxSource, dySource error of the calculated position of the source in
295// the camera
296//
297// Useful formulas can be found in TDAS 00-11 and TDAS 01-05
298//
299
300void MSourcePosfromStarPos::SourcefromStar(Double_t &f,
301 TArrayD &decStar, TArrayD &raStar,
302 Double_t &decSource, Double_t &raSource,
303 Double_t &thetaTel, Double_t &phiTel,
304 TArrayD &xStar, TArrayD &yStar,
305 TArrayD &dxStar, TArrayD &dyStar,
306 Double_t &xSource, Double_t &ySource,
307 Double_t &dxSource, Double_t &dySource)
308{
309 /*
310 *fLog << "MSourcePosfromStarPos::SourcefromStar : printout in degrees" << endl;
311 *fLog << " decStar, raStar = " << decStar[0]*kRad2Deg << ", "
312 << raStar[0]*kRad2Deg << endl;
313 *fLog << " decSource, raSource = " << decSource*kRad2Deg << ", "
314 << raSource*kRad2Deg << endl;
315 *fLog << " thetaTel, phiTel = " << thetaTel*kRad2Deg << ", "
316 << phiTel*kRad2Deg << endl;
317 *fLog << " xStar, yStar = " << xStar[0]*fMm2Deg << ", "
318 << yStar[0]*fMm2Deg << endl;
319
320 *fLog << "MSourcePosfromStarPos::SourcefromStar : printout in radians and mm" << endl;
321 *fLog << " decStar, raStar = " << decStar[0] << ", "
322 << raStar[0] << endl;
323 *fLog << " decSource, raSource = " << decSource << ", "
324 << raSource << endl;
325 *fLog << " thetaTel, phiTel = " << thetaTel << ", "
326 << phiTel << endl;
327 *fLog << " xStar, yStar = " << xStar[0] << ", "
328 << yStar[0] << endl;
329 */
330
331 // the units are assumed to be radians for theta, phi, dec and ra
332 // and mm for f, x and y
333
334
335 // calculate rotation angle alpha of sky image in camera
336 // (see TDAS 00-11, eqs. (18) and (20))
337 // a1 = cos(Lat), a3 = -sin(Lat), where Lat is the geographical latitude of La Palma
338 Double_t a1 = 0.876627;
339 Double_t a3 = -0.481171;
340
341 Double_t denom = 1./ sqrt( sin(thetaTel)*sin(phiTel) * sin(thetaTel)*sin(phiTel) +
342 ( a1*cos(thetaTel)+a3*sin(thetaTel)*cos(phiTel) ) *
343 ( a1*cos(thetaTel)+a3*sin(thetaTel)*cos(phiTel) ) );
344 Double_t cosal = - (a3 * sin(thetaTel) + a1 * cos(thetaTel) * cos(phiTel)) * denom;
345 Double_t sinal = a1 * sin(phiTel) * denom;
346
347 fObservatory->RotationAngle(thetaTel, phiTel, sinal, cosal);
348
349 //*fLog << "thetaTel, phiTel, cosal, sinal = " << thetaTel << ", "
350 // << phiTel << ", " << cosal << ", " << sinal << endl;
351
352
353 // calculate coordinates of source in system B (see TDAS 00-11, eqs. (2))
354 // note that right ascension and hour angle go into opposite directions
355 Double_t xB0 = cos(decSource) * cos(-raSource);
356 Double_t yB0 = cos(decSource) * sin(-raSource);
357 Double_t zB0 = -sin(decSource);
358
359 //*fLog << "xB0, yB0, zB0 = " << xB0 << ", " << yB0 << ", "
360 // << zB0 << endl;
361
362 //-----------------------------------------------------
363 // loop over stars
364 Double_t sumx = 0.0;
365 Double_t sumy = 0.0;
366 Double_t sumwx = 0.0;
367 Double_t sumwy = 0.0;
368
369 for (Int_t i=0; i<decStar.GetSize(); i++)
370 {
371 // calculate weights
372 Double_t weightx = 1.0 / (dxStar[i]*dxStar[i]);
373 Double_t weighty = 1.0 / (dyStar[i]*dyStar[i]);
374 sumwx += weightx;
375 sumwy += weighty;
376
377 //*fLog << "weightx, weighty = " << weightx << ", " << weighty << endl;
378
379 // calculate coordinates of star in system B (see TDAS 00-11, eqs. (2))
380 // note that right ascension and hour angle go into opposite directions
381 Double_t xB = cos(decStar[i]) * cos(-raStar[i]);
382 Double_t yB = cos(decStar[i]) * sin(-raStar[i]);
383 Double_t zB = -sin(decStar[i]);
384
385
386 //*fLog << "xB, yB, zB = " << xB << ", " << yB << ", "
387 // << zB << endl;
388
389
390 // calculate coordinates of star in a system with the basis vectors e1, e2, e3
391 // where e1 is in the direction (r0 x a)
392 // e2 is in the direction (e1 x r0)
393 // and e3 is in the direction -r0;
394 // r0 is the direction to the source
395 // and a is the earth rotation axis (pointing to the celestial north pole)
396 //
397 Double_t x = (-xB*yB0 + xB0*yB) / sqrt( xB0*xB0 + yB0*yB0 );
398 Double_t y = ( xB*xB0*zB0 + yB*yB0*zB0 - zB*(xB0*xB0 + yB0*yB0) )
399 / sqrt( xB0*xB0 + yB0*yB0 );
400 Double_t z = -(xB*xB0 + yB*yB0 + zB*zB0);
401
402 //*fLog << "x, y, z = " << x << ", " << y << ", "
403 // << z << endl;
404
405
406 // calculate coordinates of star in camera
407 Double_t xtilde = -f/z * (cosal*x - sinal*y);
408 Double_t ytilde = -f/z * (sinal*x + cosal*y);
409
410 //*fLog << "xtilde, ytile = " << xtilde << ", " << ytilde << endl;
411
412
413 // calculate coordinates of source in camera
414 // note : in real camera signs are inverted (therefore s = -1.0)
415 Double_t s = -1.0;
416 sumx += s * (s*xStar[i] - xtilde) * weightx;
417 sumy += s * (s*yStar[i] - ytilde) * weighty;
418 }
419 //-----------------------------------------------------
420
421 xSource = sumx / sumwx;
422 ySource = sumy / sumwy;
423 dxSource = 1.0 / sqrt(sumwx);
424 dySource = 1.0 / sqrt(sumwy);
425
426 /*
427 Int_t run = fRun->GetRunNumber();
428 *fLog << all << "MSourcePosfromStarPos::SourcefromStar; run, xSource, ySource = "
429 << run << " : "
430 << xSource << " +- " << dxSource << ", "
431 << ySource << " +- " << dySource << endl;
432 */
433}
434
435// --------------------------------------------------------------------------
436//
437// Get the source position and put it into MSrcPosCam
438//
439//
440Bool_t MSourcePosfromStarPos::ReInit(MParList *pList)
441{
442
443 Int_t run = fRun->GetRunNumber();
444 *fLog << all << "MSourcePosfromStarPos::ReInit; run = " << run << endl;
445
446 // define default values for source position and (theta, phi)
447 Double_t xSource_def = 51.0;
448 Double_t ySource_def = -59.0;
449 Double_t theta_def = 25.0;
450 Double_t phi_def = 90.0;
451
452 //-------------------------------------------------------------------
453 // search this run in the list
454 for (Int_t i=0; i<fSize; i++)
455 {
456 if (run == fRunNr[i])
457 {
458 //-----------------------------------------
459 // put the zenith angle into MMcEvt
460
461 Double_t thetarad = fThetaTel[i];
462 Double_t phirad = fPhiTel[i];
463 fMcEvt->SetTelescopeTheta(thetarad);
464 fMcEvt->SetTelescopePhi(phirad);
465 fMcEvt->SetReadyToSave();
466
467 *fLog << all << "theta, phi = " << thetarad*kRad2Deg << ", "
468 << phirad*kRad2Deg << " deg" << endl;
469
470 //-----------------------------------------
471 // Get source position and put it into MSrcPosCam
472
473
474 if (fStars > 0)
475 {
476 TArrayD xStar(fxStar.GetNrows());
477 TArrayD dxStar(fdxStar.GetNrows());
478 TArrayD yStar(fyStar.GetNrows());
479 TArrayD dyStar(fdyStar.GetNrows());
480 for (Int_t j=0; j<fxStar.GetNrows(); j++)
481 {
482 xStar[j] = fxStar(j, i);
483 dxStar[j] = fdxStar(j, i);
484 yStar[j] = fyStar(j, i);
485 dyStar[j] = fdyStar(j, i);
486 }
487
488 MSourcePosfromStarPos::SourcefromStar( fDistCameraReflector,
489 fDecStar, fRaStar, fDecSource, fRaSource,
490 fThetaTel[i], fPhiTel[i],
491 xStar, yStar,
492 dxStar, dyStar,
493 fxSource, fySource,
494 fdxSource, fdySource);
495
496 fSrcPos->SetXY(fxSource, fySource);
497
498 *fLog << all << "MSourcePosfromStarPos::ReInit; fRunNr, fxSource, fySource = "
499 << fRunNr[i] << ", " << fxSource << " +- " << fdxSource
500 << ", " << fySource << " +- " << fdySource << endl;
501
502 fSrcPos->SetReadyToSave();
503 }
504 else
505 {
506 // set default values
507 fxSource = xSource_def;
508 fySource = ySource_def;
509 fSrcPos->SetXY(fxSource, fySource);
510 fSrcPos->SetReadyToSave();
511
512 *fLog << warn << "MSourcePosfromStarPos::ReInit; no information on source position for run number = "
513 << run << endl;
514 *fLog << warn << " set xSource, ySource = " << fxSource << ", "
515 << fySource << " mm" << endl;
516 }
517
518
519 return kTRUE;
520 }
521 }
522 //-------------------------------------------------------------------
523
524 // set default values
525 fxSource = xSource_def;
526 fySource = ySource_def;
527 fSrcPos->SetXY(fxSource, fySource);
528 fSrcPos->SetReadyToSave();
529
530 Double_t thetadeg = theta_def;
531 Double_t thetarad = thetadeg / kRad2Deg;
532 fMcEvt->SetTelescopeTheta(thetarad);
533
534 Double_t phideg = phi_def;
535 Double_t phirad = phideg / kRad2Deg;
536 fMcEvt->SetTelescopePhi(phirad);
537 fMcEvt->SetReadyToSave();
538
539 *fLog << warn << "MSourcePosfromStarPos::ReInit; no information on theta, phi and source position for run number = "
540 << run << endl;
541 *fLog << warn << " set xSource, ySource = " << fxSource << ", "
542 << fySource << " mm" << endl;
543 *fLog << warn << " set theta, phi = " << thetadeg << ", "
544 << phideg << " deg" << endl;
545
546
547 return kTRUE;
548}
549
550// --------------------------------------------------------------------------
551//
552//
553Int_t MSourcePosfromStarPos::Process()
554{
555 //Int_t run = fRun->GetRunNumber();
556 //*fLog << "MSourcePosfromStarPos::Process; run = " << run << endl;
557
558
559 return kTRUE;
560}
561
562// --------------------------------------------------------------------------
563//
564//
565Int_t MSourcePosfromStarPos::PostProcess()
566{
567
568 return kTRUE;
569}
570
571// --------------------------------------------------------------------------
572//
573// read the data from the ASCII file and store them
574//
575void MSourcePosfromStarPos::FixSize()
576{
577 *fLog << "MSourcePosfromStarPos::FixSize; fix size of arrays : fRuns = "
578 << fRuns << ", fStars = " << fStars << endl;
579
580 fSize = fRuns;
581
582 fRunNr.Set(fSize);
583
584 fThetaTel.Set(fSize);
585 fPhiTel.Set(fSize);
586 fdThetaTel.Set(fSize);
587 fdPhiTel.Set(fSize);
588
589 Int_t fRows = fxStar.GetNrows();
590 fxStar.ResizeTo(fRows, fSize);
591 fyStar.ResizeTo(fRows, fSize);
592 fdxStar.ResizeTo(fRows, fSize);
593 fdyStar.ResizeTo(fRows, fSize);
594}
595
596// --------------------------------------------------------------------------
597//
598// read the data from the ASCII file and store them
599//
600void MSourcePosfromStarPos::ReadData()
601{
602 Float_t val;
603 Int_t ival;
604
605 // extend size of arrays if necessary
606 if ( fRuns >= (fSize-1) )
607 {
608 fSize += 100;
609
610 fRunNr.Set(fSize);
611
612 fThetaTel.Set(fSize);
613 fPhiTel.Set(fSize);
614 fdThetaTel.Set(fSize);
615 fdPhiTel.Set(fSize);
616
617 Int_t fRows = fxStar.GetNrows();
618 fxStar.ResizeTo(fRows, fSize);
619 fyStar.ResizeTo(fRows, fSize);
620 fdxStar.ResizeTo(fRows, fSize);
621 fdyStar.ResizeTo(fRows, fSize);
622 }
623
624 //-------------------
625 // read header line
626 //*fIn >> val;
627
628 //*fLog << "val =" << val << endl;
629
630 //*fIn >> val;
631 //*fIn >> val;
632 //*fIn >> val;
633 //*fIn >> val;
634 //*fIn >> val;
635 //*fIn >> val;
636 //*fIn >> val;
637 //*fIn >> val;
638 //*fIn >> val;
639 //*fIn >> val;
640 //-------------------
641
642
643 fRuns += 1;
644
645 *fIn >> ival;
646
647 //*fLog << fRuns <<"th run : " << ival << endl;
648
649 fRunNr.AddAt(ival, fRuns-1);
650
651 //*fLog << "check : fRuns, fRunNr[fRuns-1], fRunNr[fRuns] = " << fRuns << ", "
652 // << fRunNr[fRuns-1] << ", " << fRunNr[fRuns] << endl;
653
654
655 // read mjdS, hmsS, mjdE, hmsE
656 // these data are present only for ON data (fStars > 0)
657 if (fStars > 0)
658 {
659 *fIn >> val;
660 *fIn >> val;
661 *fIn >> val;
662 *fIn >> val;
663
664 *fIn >> val;
665 *fIn >> val;
666 *fIn >> val;
667 *fIn >> val;
668 }
669
670 *fIn >> val;
671 fThetaTel.AddAt(val/kRad2Deg, fRuns-1);
672 //*fLog << "val, fThetaTel[fRuns-1] = " << val << ", "
673 // << fThetaTel[fRuns-1] << endl;
674
675
676 *fIn >> val;
677 fPhiTel.AddAt(val/kRad2Deg, fRuns-1);
678 //*fLog << "val, fPhiTel[fRuns-1] = " << val << ", "
679 // << fPhiTel[fRuns-1] << endl;
680
681
682 //*fIn >> val;
683 //fdThetaTel.AddAt(val/kRad2Deg, fRuns-1);
684 //*fIn >> val;
685 //fdPhiTel.AddAt(val/kRad2Deg, fRuns-1);
686
687 // input is in [deg], convert to [mm]
688
689 //*fLog << "ReadData : fStars = " << fStars << endl;
690 for (Int_t i=0; i<fStars; i++)
691 {
692 *fIn >> val;
693 fxStar(i, fRuns-1) = val / fMm2Deg;;
694 //*fLog << "val, fxStar(i, fRuns-1) = " << val << ", "
695 // << fxStar(i, fRuns-1) << endl;
696
697 *fIn >> val;
698 fyStar(i, fRuns-1) = val / fMm2Deg;
699 //*fLog << "val, fyStar(i, fRuns-1) = " << val << ", "
700 // << fyStar(i, fRuns-1) << endl;
701
702
703 *fIn >> val;
704 //*fLog << "y=dxStar = " << val << endl;
705
706 fdxStar(i, fRuns-1) = val / fMm2Deg;
707 *fIn >> val;
708 //*fLog << "y=dyStar = " << val << endl;
709
710 fdyStar(i, fRuns-1) = val / fMm2Deg;
711 }
712
713}
714
715// --------------------------------------------------------------------------
716//
717// Add this file as the last entry in the chain
718//
719Int_t MSourcePosfromStarPos::AddFile(const char *txt, Int_t)
720{
721 TNamed *name = new TNamed(txt, "");
722 fFileNames->AddLast(name);
723
724 *fLog << "MSourcePosfromStarPos::AddFile; add file '" << txt << "'"
725 << endl;
726
727 return 1;
728}
729
730// --------------------------------------------------------------------------
731//
732// This opens the next file in the list and deletes its name from the list.
733//
734Bool_t MSourcePosfromStarPos::OpenNextFile()
735{
736 //
737 // open the input stream and check if it is really open (file exists?)
738 //
739 if (fIn)
740 delete fIn;
741 fIn = NULL;
742
743 //
744 // Check for the existence of a next file to read
745 //
746 TNamed *file = (TNamed*)fFileNames->First();
747 if (!file)
748 return kFALSE;
749
750 //
751 // open the file which is the first one in the chain
752 //
753 const char *name = file->GetName();
754
755 const char *expname = gSystem->ExpandPathName(name);
756 fIn = new ifstream(expname);
757 delete [] expname;
758
759 const Bool_t noexist = !(*fIn);
760
761
762 if (noexist)
763 *fLog << dbginf << "Cannot open file '" << name << "'" << endl;
764 else
765 *fLog << "Open file: '" << name << "'" << endl;
766
767 //
768 // Remove this file from the list of pending files
769 //
770 fFileNames->Remove(file);
771
772 return !noexist;
773}
774
775// --------------------------------------------------------------------------
776
777
778
779
780
781
782
783
784
785
786
Note: See TracBrowser for help on using the repository browser.