Index: /trunk/MagicSoft/Mars/Changelog =================================================================== --- /trunk/MagicSoft/Mars/Changelog (revision 3365) +++ /trunk/MagicSoft/Mars/Changelog (revision 3366) @@ -11,4 +11,5 @@ + 2004/03/01: Thomas Bretz @@ -20,4 +21,15 @@ - enhanced interpretation of time - fixed a bug in determin the file-type + + * mastro/MAstro.[h,cc]: + - added Coordinate2Angle + - added AngularDistance + + * mastro/MObservatory.[h,cc]: + - added data members for sin-/cos-component of latitude + - added calculation of rotation angle + + * mhist/MHStarMap.[h,cc]: + - changes to use MObservatory member function Index: /trunk/MagicSoft/Mars/mastro/MAstro.cc =================================================================== --- /trunk/MagicSoft/Mars/mastro/MAstro.cc (revision 3365) +++ /trunk/MagicSoft/Mars/mastro/MAstro.cc (revision 3366) @@ -31,4 +31,6 @@ #include "MAstro.h" +#include // TVector3 + ClassImp(MAstro); @@ -183,4 +185,12 @@ } +// -------------------------------------------------------------------------- +// +// Interpretes a string ' - 12 30 00.0' or '+ 12 30 00.0' +// as floating point value -12.5 or 12.5. If interpretation is +// successfull kTRUE is returned, otherwise kFALSE. ret is not +// touched if interpretation was not successfull. The successfull +// interpreted part is removed from the TString. +// Bool_t MAstro::String2Angle(TString &str, Double_t &ret) { @@ -202,4 +212,31 @@ } +// -------------------------------------------------------------------------- +// +// Interpretes a string '-12:30:00.0', '12:30:00.0' or '+12:30:00.0' +// as floating point value -12.5, 12.5 or 12.5. If interpretation is +// successfull kTRUE is returned, otherwise kFALSE. ret is not +// touched if interpretation was not successfull. +// +Bool_t MAstro::Coordinate2Angle(const TString &str, Double_t &ret) +{ + Char_t sgn = str[0]=='-' ? '-' : '+'; + Int_t d; + UInt_t m; + Float_t s; + + const int n=sscanf(str[0]=='+'||str[0]=='-' ? str.Data()+1 : str.Data(), "%d:%d:%f", &d, &m, &s); + + if (n!=3) + return kFALSE; + + ret = Dms2Deg(d, m, s, sgn); + return kTRUE; +} + +// -------------------------------------------------------------------------- +// +// Return year y, month m and day d corresponding to Mjd. +// void MAstro::Mjd2Ymd(UInt_t mjd, UShort_t &y, Byte_t &m, Byte_t &d) { @@ -214,4 +251,8 @@ } +// -------------------------------------------------------------------------- +// +// Return Mjd corresponding to year y, month m and day d. +// Int_t MAstro::Ymd2Mjd(UShort_t y, Byte_t m, Byte_t d) { @@ -237,2 +278,19 @@ return 1461L*lm10/4 + (306*((m+9)%12)+5)/10 - (3*((lm10+188)/100))/4 + d - 2399904; } + +// -------------------------------------------------------------------------- +// +// theta0, phi0 [rad]: polar angle/zenith distance, azimuth of 1st object +// theta1, phi1 [rad]: polar angle/zenith distance, azimuth of 2nd object +// AngularDistance [rad]: Angular distance between two objects +// +Double_t MAstro::AngularDistance(Double_t theta0, Double_t phi0, Double_t theta1, Double_t phi1) +{ + TVector3 v0(1); + v0.Rotate(phi0, theta0); + + TVector3 v1(1); + v1.Rotate(phi1, theta1); + + return v0.Angle(v1); +} Index: /trunk/MagicSoft/Mars/mastro/MAstro.h =================================================================== --- /trunk/MagicSoft/Mars/mastro/MAstro.h (revision 3365) +++ /trunk/MagicSoft/Mars/mastro/MAstro.h (revision 3366) @@ -38,5 +38,8 @@ static Bool_t String2Angle(TString &str, Double_t &ret); + static Bool_t Coordinate2Angle(const TString &str, Double_t &ret); + static Double_t AngularDistance(Double_t theta0, Double_t phi0, Double_t theta1, Double_t phi1); + static void Mjd2Ymd(UInt_t mjd, UShort_t &y, Byte_t &m, Byte_t &d); static Int_t Ymd2Mjd(UShort_t y, Byte_t m, Byte_t d); Index: /trunk/MagicSoft/Mars/mastro/MObservatory.cc =================================================================== --- /trunk/MagicSoft/Mars/mastro/MObservatory.cc (revision 3365) +++ /trunk/MagicSoft/Mars/mastro/MObservatory.cc (revision 3366) @@ -81,5 +81,5 @@ fHeight = 2196.5; // m fObservatoryName = "Observatorio del Roque de los Muchachos (Magic1)"; - return; + break; case kWuerzburgCity: @@ -88,7 +88,9 @@ fHeight = 300; fObservatoryName = "Wuerzburg City"; - return; + break; + } - } + fSinLatitude = TMath::Sin(fLatitude); + fCosLatitude = TMath::Cos(fLatitude); } @@ -102,2 +104,67 @@ } +// -------------------------------------------------------------------------- +// +// RotationAngle +// +// calculates the angle for the rotation of the sky image in the camera; +// this angle is a function of the local coordinates +// +// theta [rad]: polar angle/zenith distance +// phi [rad]: rotation angle/azimuth +// +// Return sin/cos component of angle +// +// calculate rotation angle alpha of sky image in camera +// (see TDAS 00-11, eqs. (18) and (20)) +// +void MObservatory::RotationAngle(Double_t theta, Double_t phi, Double_t &sin, Double_t &cos) const +{ + const Double_t sint = TMath::Sin(theta); + const Double_t cost = TMath::Cos(theta); + + const Double_t sinl = fSinLatitude*sint; + const Double_t cosl = fCosLatitude*cost; + + const Double_t sinp = TMath::Sin(phi); + const Double_t cosp = TMath::Cos(phi); + + const Double_t v1 = sint*sinp; + const Double_t v2 = cosl - sinl*cosp; + + const Double_t denom = TMath::Sqrt(v1*v1 + v2*v2); + + sin = (fCosLatitude*sinp) / denom; + cos = sinl + cosl*cosp / denom; +} + +// -------------------------------------------------------------------------- +// +// RotationAngle +// +// calculates the angle for the rotation of the sky image in the camera; +// this angle is a function of the local coordinates +// +// theta [rad]: polar angle/zenith distance +// phi [rad]: rotation angle/azimuth +// +// Return RotationAngle in rad +// +// calculate rotation angle alpha of sky image in camera +// (see TDAS 00-11, eqs. (18) and (20)) +// +Double_t MObservatory::RotationAngle(Double_t theta, Double_t phi) const +{ + const Double_t sint = TMath::Sin(theta); + const Double_t cost = TMath::Cos(theta); + + const Double_t sinp = TMath::Sin(phi); + const Double_t cosp = TMath::Cos(phi); + + const Double_t v1 = sint*sinp; + const Double_t v2 = fCosLatitude*cost - fSinLatitude*sint*cosp; + + const Double_t denom = TMath::Sqrt(v1*v1 + v2*v2); + + return TMath::ASin((fCosLatitude*sinp) / denom); +} Index: /trunk/MagicSoft/Mars/mastro/MObservatory.h =================================================================== --- /trunk/MagicSoft/Mars/mastro/MObservatory.h (revision 3365) +++ /trunk/MagicSoft/Mars/mastro/MObservatory.h (revision 3366) @@ -23,4 +23,7 @@ Double_t fLatitude; //! [rad] Latitude of observatory (+ north) + Double_t fSinLatitude; //! Sin component for faster access + Double_t fCosLatitude; //! Cos component for faster access + Double_t fHeight; //! [m] height of observatory @@ -46,5 +49,11 @@ Double_t GetElong() const { return fLongitude; } //[rad] + Double_t GetSinPhi() const { return fSinLatitude; } + Double_t GetCosPhi() const { return fCosLatitude; } + Double_t GetHeight() const { return fHeight; } + + void RotationAngle(Double_t theta, Double_t phi, Double_t &sin, Double_t &cos) const; + Double_t RotationAngle(Double_t theta, Double_t phi) const; LocationName_t GetObservatoryKey() const { return fObservatoryKey; } Index: /trunk/MagicSoft/Mars/mhist/MHStarMap.cc =================================================================== --- /trunk/MagicSoft/Mars/mhist/MHStarMap.cc (revision 3365) +++ /trunk/MagicSoft/Mars/mhist/MHStarMap.cc (revision 3366) @@ -17,5 +17,5 @@ ! ! Author(s): Thomas Bretz 12/2000 -! Wolfgang Wittek 02/2004 +! Author(s): Wolfgang Wittek 02/2004 ! ! Copyright: MAGIC Software Development, 2000-2004 @@ -178,6 +178,6 @@ Bool_t MHStarMap::ReInit(MParList *pList) { - MObservatory *obs = (MObservatory*)pList->FindObject(AddSerialNumber("MObservatory")); - if (!obs) + fObservatory = (MObservatory*)pList->FindObject(AddSerialNumber("MObservatory")); + if (!fObservatory) { *fLog << err << "MObservatory not found... aborting" << endl; @@ -188,7 +188,4 @@ fCosLat = TMath::Cos(obs->GetPhi()); - *fLog << "fSinLat, fCosLat = " << fSinLat << ", " - << fCosLat << endl; - return kTRUE; } @@ -196,29 +193,10 @@ // -------------------------------------------------------------------------- // -// RotationAngle -// -// calculates the angle for the rotation of the sky image in the camera; -// this angle is a function of the local coordinates (theta, phi) in rad. -// -// calculate rotation angle alpha of sky image in camera -// (see TDAS 00-11, eqs. (18) and (20)) -// -void MHStarMap::GetRotationAngle(const Double_t &theta, const Double_t &phi, - Double_t &cosal, Double_t &sinal) -{ - - const Double_t sint = TMath::Sin(theta); - const Double_t cost = TMath::Cos(theta); - - const Double_t sinp = TMath::Sin(phi); - const Double_t cosp = TMath::Cos(phi); - - const Double_t v1 = sint*sinp; - const Double_t v2 = fCosLat*cost - fSinLat*sint*cosp; - - const Double_t denom = 1./ TMath::Sqrt(v1*v1 + v2*v2); - - cosal = (fSinLat * sint - fCosLat * cost * cosp) * denom; - sinal = fCosLat * sinp * denom; +// Calls MObservatory::RotationAngle +// +void MHStarMap::GetRotationAngle(Double_t &sin, Double_t &cos) +{ + fObservatory->RotationAngle(fMcEvt->GetTelescopeTheta(), + fMcEvt->GetTelescopePhi(), sin, cos); } @@ -270,5 +248,5 @@ Double_t cosal; Double_t sinal; - GetRotationAngle(theta, phi, cosal, sinal); + GetRotationAngle(sinal, cosal); if (m>-1 && m<1) Index: /trunk/MagicSoft/Mars/mhist/MHStarMap.h =================================================================== --- /trunk/MagicSoft/Mars/mhist/MHStarMap.h (revision 3365) +++ /trunk/MagicSoft/Mars/mhist/MHStarMap.h (revision 3366) @@ -10,10 +10,12 @@ class MSrcPosCam; class MMcEvt; +class MObservatory; class MHStarMap : public MH { private: - MSrcPosCam *fSrcPos; //! - MMcEvt *fMcEvt; //! + MSrcPosCam *fSrcPos; //! + MMcEvt *fMcEvt; //! + MObservatory *fObservatory; //! TH2F *fStarMap; //-> @@ -23,11 +25,9 @@ Bool_t fUseMmScale; - Float_t fCosLat; //! - Float_t fSinLat; //! - void PrepareDrawing() const; void Paint(Option_t *opt=""); + void GetRotationAngle(Double_t &sinangle, Double_t &cosangle); Bool_t SetupFill(const MParList *pList); @@ -46,7 +46,4 @@ TH2F *GetHist() { return fStarMap; } - void GetRotationAngle(const Double_t &thetatel, const Double_t &phitel, - Double_t &cosangle, Double_t &sinangle); - void Draw(Option_t *opt=NULL); TObject *DrawClone(Option_t *opt=NULL) const;