Main Page | Modules | Namespace List | Class Hierarchy | Alphabetical List | Class List | Directories | File List | Namespace Members | Class Members | File Members
AstroUtil

AstCoordinate.h File Reference

#include "Conventions/DetectorType.h"

Go to the source code of this file.

Namespaces

namespace  AstUtil

Functions

void LocalToIdeal (double dcosx, double dcosy, double dcosz, DetectorType::Detector_t dettype, double &dcosx_ideal, double &dcosy_ideal, double &dcosz_ideal)
void IdealToLocal (double dcosx_ideal, double dcosy_ideal, double dcosz_ideal, DetectorType::Detector_t dettype, double &dcosx, double &dcosy, double &dcosz)
void HorizonToLocal (double altitude, double azimuth, DetectorType::Detector_t dettype, double &dcosx, double &dcosy, double &dcosz)
void LocalToHorizon (double dcosx, double dcosy, double dcosz, DetectorType::Detector_t dettype, double &altitude, double &azimuth)
void HorizonToIdeal (double altitude, double azimuth, double &dcosx_ideal, double &dcosy_ideal, double &dcosz_ideal)
void IdealToHorizon (double dcosx_ideal, double dcosy_ideal, double dcosz_ideal, double &altitude, double &azimuth)
void HorizonToEquatorial (double altitude, double azimuth, double latitude, double &hourangle, double &declination)
void EquatorialToHorizon (double hourangle, double declination, double latitude, double &altitude, double &azimuth)
void CelestialToEquatorial (double ra, double gmst, double longitude, double &hourangle)
void EquatorialToCelestial (double hourangle, double gmst, double longitude, double &ra)
void CelestialToEcliptic (double declination, double ra, double &beta, double &lamda)
void EclipticToCelestial (double beta, double lamda, double &declination, double &ra)

Function Documentation

void AstUtil::CelestialToEcliptic (  double  declination,
double  ra,
double &  beta,
double &  lamda
) 
 

Definition at line 329 of file AstCoordinate.cxx.

00331 { // Static method to convert from celestial to ecliptic coordinates
00332 
00333     const double epsilon  = 23.43;
00334     
00335     double sindelta = sin(declination);
00336     double sinalpha = sin(ra);
00337     double cosdelta = cos(declination);
00338     double cosalpha = cos(ra);
00339     double cosbeta  = cos(beta);
00340     
00341     //We will need both trigonometric numbers of epsilon, the constant angle between
00342     //the ecliptic and the celestial equator great circles. 
00343     const double Alpha = cos(Mphysical::pi*epsilon/180.);
00344     const double Beta  = sin(Mphysical::pi*epsilon/180.);
00345     
00346     beta  =  asin(sindelta*Alpha-cosdelta*Beta*sinalpha);
00347     if (cosbeta!=0){
00348         lamda =  acos(cosalpha*cosdelta/cosbeta); 
00349     } else { //the point lies on pole of ecliptic, and lamda has no meaning
00350         lamda = 0;
00351     }
00352     if (lamda < 0.) lamda += 360.;
00353     return;
00354 }

void AstUtil::CelestialToEquatorial (  double  ra,
double  gmst,
double  longitude,
double &  hourangle
) 
 

Definition at line 302 of file AstCoordinate.cxx.

References AstUtil::GSTToLST().

Referenced by AstValidate::TestEquatorialToCelestial().

00303                                                               {
00304   // Static method to convert from celestial to equatorial coordinates
00305 
00306   double lst;
00307   AstUtil::GSTToLST(gst,longitude,lst);
00308   lst *= 15.; // convert from hours to degrees
00309   hourangle = lst - ra;
00310   if ( hourangle < 0. ) hourangle += 360.;
00311  
00312   return;
00313 
00314 }

void AstUtil::EclipticToCelestial (  double  beta,
double  lamda,
double &  declination,
double &  ra
) 
 

Definition at line 357 of file AstCoordinate.cxx.

00359 { // Static method to convert from ecliptic to celestial coordinates
00360 
00361     const double epsilon  = 23.43;
00362     
00363     double sinlamda = sin(lamda);
00364     double coslamda = cos(lamda);
00365     double cosbeta  = cos(beta);
00366     double sinbeta  = sin(beta);
00367     double cosdelta = cos(declination);
00368     
00369     //We will need both trigonometric numbers of epsilon, the constant angle between
00370     //the ecliptic and the celestial equator great circles. 
00371     const double Alpha = cos(Mphysical::pi*epsilon/180.);
00372     const double Beta  = sin(Mphysical::pi*epsilon/180.);
00373     
00374     declination  =  asin(sinbeta*Alpha+cosbeta*Beta*sinlamda);
00375     if (cosdelta!=0){
00376         ra =  acos(coslamda*cosbeta/cosdelta); 
00377     } else { //the point lies on pole of equator, and ra has no meaning
00378         ra = 0;
00379     }
00380     if (ra <0.) ra += 360.;
00381     return;
00382 }

void AstUtil::EquatorialToCelestial (  double  hourangle,
double  gmst,
double  longitude,
double &  ra
) 
 

Definition at line 316 of file AstCoordinate.cxx.

References AstUtil::GSTToLST().

Referenced by NtpSRModule::FillNtpTrackCosmicRay(), and AstValidate::TestEquatorialToCelestial().

00317                                                                   {
00318   // Static method to convert from equatorial to celestial coordinates
00319 
00320   double lst; 
00321   AstUtil::GSTToLST(gst,longitude,lst);
00322   lst *= 15.; // convert from hours to degrees
00323   ra = lst - hourangle; // in degrees
00324   if ( ra < 0. ) ra += 360.;
00325 
00326   return;
00327 }

void AstUtil::EquatorialToHorizon (  double  hourangle,
double  declination,
double  latitude,
double &  altitude,
double &  azimuth
) 
 

Definition at line 243 of file AstCoordinate.cxx.

Referenced by AstValidate::TestHorizonToEquatorial().

00244                                                         {
00245   // Static method to convert from equatorial to horizon coordinates
00246   // From Practical Ephemeris Calculations by O. Montenbruck:
00247   // 1)sin(alt)          =  sin(lat)*sin(dec) + cos(lat)*cos(dec)*cos(ha)
00248   // 2)cos(alt)*sin(azi) = -cos(dec)*sin(ha)
00249   // 3)cos(alt)*cos(azi) =  cos(lat)*sin(dec) - sin(lat)*cos(dec)*cos(ha)
00250   // This implementation is adapted from an M. Thomson Soudan2 fortran
00251   // routine.
00252 
00253   double latrad = latitude*Mphysical::pi/180.;
00254   double harad  = hourangle*Mphysical::pi/180.;
00255   double decrad = declination*Mphysical::pi/180.;
00256 
00257   double sinlat = sin(latrad);
00258   double sinha  = sin(harad);
00259   double sindec = sin(decrad);
00260   double coslat = cos(latrad);
00261   double cosha  = cos(harad);
00262   double cosdec = cos(decrad);
00263 
00264   // Solve equation 1 for altitude
00265   double sinalt = sindec*sinlat + cosdec*cosha*coslat;
00266   if ( sinalt >= 1 ) altitude = Mphysical::pi/2.;
00267   else if ( sinalt <= -1 ) altitude = -Mphysical::pi/2.;
00268   else altitude = asin(sinalt);
00269 
00270   // Solve equation 2 to determine azimuth in range -pi/2 to pi/2.
00271   double du = 90.0001*Mphysical::pi/180.;
00272   double dl = 89.9999*Mphysical::pi/180.;
00273   if (fabs(altitude) > dl && fabs(altitude) < du ) {
00274     azimuth = 0;
00275   }
00276   else {
00277     double sinazi = -cosdec*sinha/cos(altitude);
00278     if ( sinazi >= 1 ) azimuth = Mphysical::pi/2.;
00279     else if (sinazi <= -1 ) azimuth = -Mphysical::pi/2.; 
00280     else azimuth = asin(sinazi);
00281   }
00282 
00283   // Use 3rd equation to fully determine azimuth
00284   double lhs = cos(altitude)*cos(azimuth);
00285   double rhs = sindec*coslat - cosdec*cosha*sinlat;
00286 
00287   double prod = lhs*rhs;
00288   if ( prod != fabs(prod) ) {
00289     if ( azimuth > 0. ) azimuth = Mphysical::pi - azimuth;
00290     else azimuth = -Mphysical::pi - azimuth;
00291   }
00292 
00293   if ( azimuth < 0. ) azimuth += 2.*Mphysical::pi;
00294 
00295   altitude = altitude*180./Mphysical::pi;  // degrees
00296   azimuth  = azimuth*180./Mphysical::pi;   // degrees
00297 
00298   return;
00299 
00300 }

void AstUtil::HorizonToEquatorial (  double  altitude,
double  azimuth,
double  latitude,
double &  hourangle,
double &  declination
) 
 

Definition at line 184 of file AstCoordinate.cxx.

Referenced by NtpSRModule::FillNtpTrackCosmicRay(), and AstValidate::TestHorizonToEquatorial().

00185                                                                        {
00186   // Static method to convert from horizon to equatorial coordinates.
00187   // From Practical Ephemeris Calculations by O. Montenbruck:
00188   // 1)sin(dec)         = sin(lat)*sin(alt) + cos(lat)*cos(alt)*cos(azi)
00189   // 2)cos(dec)*sin(ha) = -cos(alt)*sin(azi)
00190   // 3)cos(dec)*cos(ha) = cos(lat)*sin(alt)-sin(lat)*cos(alt)*cos(azi)
00191   // This implementation is adapted from an M. Thomson Soudan2 fortran
00192   // routine.
00193 
00194   double latrad = latitude*Mphysical::pi/180.;
00195   double altrad = altitude*Mphysical::pi/180.;
00196   double azirad = azimuth*Mphysical::pi/180.;
00197 
00198   double sinlat = sin(latrad);
00199   double sinalt = sin(altrad);
00200   double sinazi = sin(azirad);
00201   double coslat = cos(latrad);
00202   double cosalt = cos(altrad);
00203   double cosazi = cos(azirad);
00204 
00205   // solve equation 1 to determine declination
00206   double sindec = sinalt*sinlat + cosalt*cosazi*coslat;
00207   if ( sindec >= 1 ) declination = Mphysical::pi/2.;
00208   else if ( sindec <= -1 ) declination = -Mphysical::pi/2.;
00209   else declination = asin(sindec);
00210 
00211   // solve equation 2 to determine hourangle in range -pi/2 to pi/2
00212   double du = 90.0001*Mphysical::pi/180.;
00213   double dl = 89.9999*Mphysical::pi/180.;
00214   if ( fabs(declination) > dl && fabs(declination) < du ) {
00215     hourangle = 0.;
00216   }
00217   else {
00218    double sinhourangle = -cosalt*sinazi/cos(declination);
00219    if ( sinhourangle >= 1 ) hourangle = Mphysical::pi/2.;
00220    else if ( sinhourangle <= -1 ) hourangle = -Mphysical::pi/2.;
00221    else hourangle = asin(sinhourangle);
00222   }
00223 
00224   // use 3rd equation to fully determine hourangle
00225   double lhs = cos(declination)*cos(hourangle);
00226   double rhs = sinalt*coslat - cosalt*cosazi*sinlat;
00227 
00228   double prod = lhs*rhs;
00229   if ( prod != fabs(prod) ) {
00230     if ( hourangle > 0. ) hourangle = Mphysical::pi - hourangle;
00231     else hourangle = -Mphysical::pi - hourangle;
00232   }
00233 
00234   if ( hourangle < 0. ) hourangle += 2.*Mphysical::pi;
00235 
00236   hourangle   = hourangle*180./Mphysical::pi; // degrees
00237   declination = declination*180./Mphysical::pi; // degrees
00238 
00239   return;
00240 
00241 }

void AstUtil::HorizonToIdeal (  double  altitude,
double  azimuth,
double &  dcosx_ideal,
double &  dcosy_ideal,
double &  dcosz_ideal
) 
 

Definition at line 158 of file AstCoordinate.cxx.

Referenced by AstUtil::HorizonToLocal(), and AstValidate::TestIdealToHorizon().

00159                                                                    {
00160   // Static method to convert horizon coordinates to ideal detector
00161   // coordinates.
00162   double zenrad = (90. - altitude)*Mphysical::pi/180.; 
00163   double azirad = azimuth*Mphysical::pi/180.;
00164   dcosy_ideal =  cos(zenrad);
00165   dcosx_ideal = -sin(zenrad)*sin(azirad);
00166   dcosz_ideal =  sin(zenrad)*cos(azirad);
00167   return;
00168 
00169 }

void AstUtil::HorizonToLocal (  double  altitude,
double  azimuth,
DetectorType::Detector_t  dettype,
double &  dcosx,
double &  dcosy,
double &  dcosz
) 
 

Definition at line 171 of file AstCoordinate.cxx.

References AstUtil::HorizonToIdeal(), and AstUtil::IdealToLocal().

00173                          {
00174   // Static method to convert horizon coordinates to direction cosines
00175 
00176   double dcosx_ideal,dcosy_ideal,dcosz_ideal;
00177   AstUtil::HorizonToIdeal(altitude,azimuth,dcosx_ideal,dcosy_ideal,
00178                                                        dcosz_ideal);
00179   AstUtil::IdealToLocal(dcosx_ideal,dcosy_ideal,dcosz_ideal,dettype,
00180                         dcosx,dcosy,dcosz);
00181   return;
00182 }

void AstUtil::IdealToHorizon (  double  dcosx_ideal,
double  dcosy_ideal,
double  dcosz_ideal,
double &  altitude,
double &  azimuth
) 
 

Definition at line 146 of file AstCoordinate.cxx.

Referenced by AstUtil::LocalToHorizon(), and AstValidate::TestIdealToHorizon().

00147                                                                       {
00148   // Static method to convert ideal coordinate system direction cosines to 
00149   // horizon coordinates.
00150 
00151   double zenith = acos(dcosy_ideal)*180./Mphysical::pi;
00152   altitude = 90. - zenith;
00153   azimuth = 180.*atan2(-dcosx_ideal,dcosz_ideal)/Mphysical::pi;
00154   if ( azimuth < 0. ) azimuth += 360.;
00155 
00156 }

void AstUtil::IdealToLocal (  double  dcosx_ideal,
double  dcosy_ideal,
double  dcosz_ideal,
DetectorType::Detector_t  dettype,
double &  dcosx,
double &  dcosy,
double &  dcosz
) 
 

Definition at line 100 of file AstCoordinate.cxx.

References AstUtil::GetDetRotMatrixIdealToLocal(), and MSG.

Referenced by AstUtil::HorizonToLocal(), and AstValidate::TestLocalToIdeal().

00102                                                    {
00103   // Static method to convert direction cosines from ideal detector 
00104   // coordinates to local detector coordinates.
00105 
00106   const double* rotmatrix = AstUtil::GetDetRotMatrixIdealToLocal(dettype);
00107   if ( !rotmatrix ) {
00108     MSG("Ast",Msg::kWarning) 
00109       << "No rotation matrix available for requested det type " << dettype
00110       << ".\n No rotation from ideal to local coordinates  will be applied." 
00111       << endl;
00112     dcosx = dcosx_ideal;
00113     dcosy = dcosy_ideal;
00114     dcosz = dcosz_ideal;
00115     return;
00116   }
00117 
00118   dcosx = rotmatrix[0]*dcosx_ideal
00119         + rotmatrix[1]*dcosy_ideal
00120         + rotmatrix[2]*dcosz_ideal;
00121   dcosy = rotmatrix[3]*dcosx_ideal
00122         + rotmatrix[4]*dcosy_ideal
00123         + rotmatrix[5]*dcosz_ideal;
00124   dcosz = rotmatrix[6]*dcosx_ideal
00125         + rotmatrix[7]*dcosy_ideal
00126         + rotmatrix[8]*dcosz_ideal;
00127   return;
00128 
00129 }

void AstUtil::LocalToHorizon (  double  dcosx,
double  dcosy,
double  dcosz,
DetectorType::Detector_t  dettype,
double &  altitude,
double &  azimuth
) 
 

Definition at line 131 of file AstCoordinate.cxx.

References AstUtil::IdealToHorizon(), and AstUtil::LocalToIdeal().

Referenced by MakeAlignmentModule::Ana(), and NtpSRModule::FillNtpTrackCosmicRay().

00132                                                                        {
00133   // Static method to convert direction cosines to horizon coordinates
00134   // The direction cosines are local to the detector type specified by dettype.
00135 
00136 
00137   double dcosx_ideal,dcosy_ideal,dcosz_ideal;
00138 
00139   AstUtil::LocalToIdeal(dcosx,dcosy,dcosz,dettype,dcosx_ideal,
00140                          dcosy_ideal,dcosz_ideal);
00141   AstUtil::IdealToHorizon(dcosx_ideal,dcosy_ideal,dcosz_ideal,altitude,
00142                              azimuth);
00143 
00144 }

void AstUtil::LocalToIdeal (  double  dcosx,
double  dcosy,
double  dcosz,
DetectorType::Detector_t  dettype,
double &  dcosx_ideal,
double &  dcosy_ideal,
double &  dcosz_ideal
) 
 

Definition at line 69 of file AstCoordinate.cxx.

References AstUtil::GetDetRotMatrixLocalToIdeal(), and MSG.

Referenced by AstUtil::LocalToHorizon(), and AstValidate::TestLocalToIdeal().

00071                                                 {
00072   // Static method to convert direction cosines from local detector
00073   // coordinates to ideal detector coordinates.
00074 
00075   const double* rotmatrix = AstUtil::GetDetRotMatrixLocalToIdeal(dettype);
00076   if ( !rotmatrix ) {
00077     MSG("Ast",Msg::kWarning) 
00078       << "No rotation matrix available for requested det type " << dettype
00079       << ".\n No rotation from local to ideal coordinates will be applied." 
00080       << endl;
00081     dcosx_ideal = dcosx;
00082     dcosy_ideal = dcosy;
00083     dcosz_ideal = dcosz;
00084     return;
00085   }
00086 
00087   dcosx_ideal = rotmatrix[0]*dcosx
00088               + rotmatrix[1]*dcosy
00089               + rotmatrix[2]*dcosz;
00090   dcosy_ideal = rotmatrix[3]*dcosx
00091               + rotmatrix[4]*dcosy
00092               + rotmatrix[5]*dcosz;
00093   dcosz_ideal = rotmatrix[6]*dcosx
00094               + rotmatrix[7]*dcosy
00095               + rotmatrix[8]*dcosz;
00096   return;
00097 
00098 }

Generated on Thu Nov 1 15:54:09 2007 for loon by  doxygen 1.3.9.1