StelUtils.cpp   StelUtils.cpp 
skipping to change at line 26 skipping to change at line 26
* along with this program; if not, write to the Free Software * along with this program; if not, write to the Free Software
* Foundation, Inc., 51 Franklin Street, Suite 500, Boston, MA 02110-1335, USA. * Foundation, Inc., 51 Franklin Street, Suite 500, Boston, MA 02110-1335, USA.
*/ */
#ifdef CYGWIN #ifdef CYGWIN
#include <malloc.h> #include <malloc.h>
#endif #endif
#include "StelUtils.hpp" #include "StelUtils.hpp"
#include "VecMath.hpp" #include "VecMath.hpp"
#include <QBuffer>
#include <QString> #include <QString>
#include <QStringList> #include <QStringList>
#include <QTextStream> #include <QTextStream>
#include <QFile> #include <QFile>
#include <QDebug> #include <QDebug>
#include <QLocale> #include <QLocale>
#include <QRegExp> #include <QRegExp>
#include <QProcess> #include <QProcess>
#include <QSysInfo> #include <QSysInfo>
#include <cmath> // std::fmod #include <cmath> // std::fmod
skipping to change at line 476 skipping to change at line 477
return Vec3f(0.f,0.f,0.f); return Vec3f(0.f,0.f,0.f);
return Vec3f(s[0].toFloat(),s[1].toFloat(),s[2].toFloat()); return Vec3f(s[0].toFloat(),s[1].toFloat(),s[2].toFloat());
} }
Vec3f strToVec3f(const QString& s) Vec3f strToVec3f(const QString& s)
{ {
return strToVec3f(s.split(",")); return strToVec3f(s.split(","));
} }
Vec4d strToVec4d(const QStringList &s)
{
if(s.size()<4)
return Vec4d(0.0,0.0,0.0,0.0);
return Vec4d(s[0].toDouble(), s[1].toDouble(), s[2].toDouble(), s[3]
.toDouble());
}
Vec4d strToVec4d(const QString& str)
{
return strToVec4d(str.split(","));
}
QString vec3fToStr(const Vec3f &v)
{
return QString("%1,%2,%3")
.arg(v[0],0,'f',6)
.arg(v[1],0,'f',6)
.arg(v[2],0,'f',6);
}
QString vec4dToStr(const Vec4d &v)
{
return QString("%1,%2,%3,%4")
.arg(v[0],0,'f',10)
.arg(v[1],0,'f',10)
.arg(v[2],0,'f',10)
.arg(v[3],0,'f',10);
}
// Converts a Vec3f to HTML color notation. // Converts a Vec3f to HTML color notation.
QString vec3fToHtmlColor(const Vec3f& v) QString vec3fToHtmlColor(const Vec3f& v)
{ {
return QString("#%1%2%3") return QString("#%1%2%3")
.arg(qMin(255, int(v[0] * 255)), 2, 16, QChar('0')) .arg(qMin(255, int(v[0] * 255)), 2, 16, QChar('0'))
.arg(qMin(255, int(v[1] * 255)), 2, 16, QChar('0')) .arg(qMin(255, int(v[1] * 255)), 2, 16, QChar('0'))
.arg(qMin(255, int(v[2] * 255)), 2, 16, QChar('0')); .arg(qMin(255, int(v[2] * 255)), 2, 16, QChar('0'));
} }
Vec3f htmlColorToVec3f(const QString& c) Vec3f htmlColorToVec3f(const QString& c)
skipping to change at line 548 skipping to change at line 579
*lng = atan2(v[1],v[0]); *lng = atan2(v[1],v[0]);
} }
void rectToSphe(double *lng, double *lat, const Vec3f& v) void rectToSphe(double *lng, double *lat, const Vec3f& v)
{ {
double r = v.length(); double r = v.length();
*lat = asin(v[2]/r); *lat = asin(v[2]/r);
*lng = atan2(v[1],v[0]); *lng = atan2(v[1],v[0]);
} }
void ctRadec2Ecl(const double raRad, const double decRad, const double eclR ad, double *lambdaRad, double *betaRad) void equToEcl(const double raRad, const double decRad, const double eclRad, double *lambdaRad, double *betaRad)
{ {
*lambdaRad=std::atan2(std::sin(raRad)*std::cos(eclRad)+std::tan(decR ad)*std::sin(eclRad), std::cos(raRad)); *lambdaRad=std::atan2(std::sin(raRad)*std::cos(eclRad)+std::tan(decR ad)*std::sin(eclRad), std::cos(raRad));
*betaRad=std::asin(std::sin(decRad)*std::cos(eclRad)-std::cos(decRad )*std::sin(eclRad)*std::sin(raRad)); *betaRad=std::asin(std::sin(decRad)*std::cos(eclRad)-std::cos(decRad )*std::sin(eclRad)*std::sin(raRad));
} }
double getDecAngle(const QString& str) double getDecAngle(const QString& str)
{ {
QRegExp re1("^\\s*([\\+\\-])?\\s*(\\d+)\\s*([hHDd\xBA])\\s*(\\d+)\\s *['Mm]\\s*(\\d+(\\.\\d+)?)\\s*[\"Ss]\\s*([NSEWnsew])?\\s*$"); // DMS/HMS QRegExp re1("^\\s*([\\+\\-])?\\s*(\\d+)\\s*([hHDd\xBA])\\s*(\\d+)\\s *['Mm]\\s*(\\d+(\\.\\d+)?)\\s*[\"Ss]\\s*([NSEWnsew])?\\s*$"); // DMS/HMS
QRegExp re2("^\\s*([\\+\\-])?\\s*(\\d+(\\.\\d+)?).?([NSEWnsew])?\\s* $"); // Decimal QRegExp re2("^\\s*([\\+\\-])?\\s*(\\d+(\\.\\d+)?).?([NSEWnsew])?\\s* $"); // Decimal
QRegExp re3("([+-]?[\\d.]+)°(?:([\\d.]+)')?(?:([\\d.]+)\")?"); // DM S like +121°33'38.28" QRegExp re3("([+-]?[\\d.]+)°(?:([\\d.]+)')?(?:([\\d.]+)\")?"); // DM S like +121°33'38.28"
skipping to change at line 740 skipping to change at line 771
.arg(minute,2,10,QLatin1Char('0')) .arg(minute,2,10,QLatin1Char('0'))
.arg(second,2,10,QLatin1Char('0')); .arg(second,2,10,QLatin1Char('0'));
if (year < 0) if (year < 0)
{ {
res.prepend("-"); res.prepend("-");
} }
return res; return res;
} }
// Format the date per the fmt. // Format the date per the fmt.
QString localeDateString(const int year, const int month, const int day, co nst int dayOfWeek, const QString fmt) QString localeDateString(const int year, const int month, const int day, co nst int dayOfWeek, const QString &fmt)
{ {
/* we have to handle the year zero, and the years before qdatetime c an represent. */ /* we have to handle the year zero, and the years before qdatetime c an represent. */
const QLatin1Char quote('\''); const QLatin1Char quote('\'');
QString out; QString out;
int quotestartedat = -1; int quotestartedat = -1;
for (int i = 0; i < (int)fmt.length(); i++) for (int i = 0; i < (int)fmt.length(); i++)
{ {
if (fmt.at(i) == quote) if (fmt.at(i) == quote)
{ {
skipping to change at line 1273 skipping to change at line 1304
double getDeltaTByEspenakMeeus(const double jDay) double getDeltaTByEspenakMeeus(const double jDay)
{ {
int year, month, day; int year, month, day;
getDateFromJulianDay(jDay, &year, &month, &day); getDateFromJulianDay(jDay, &year, &month, &day);
// Note: the method here is adapted from // Note: the method here is adapted from
// "Five Millennium Canon of Solar Eclipses" [Espenak and Meeus, 200 6] // "Five Millennium Canon of Solar Eclipses" [Espenak and Meeus, 200 6]
// A summary is described here: // A summary is described here:
// http://eclipse.gsfc.nasa.gov/SEhelp/deltatpoly2004.html // http://eclipse.gsfc.nasa.gov/SEhelp/deltatpoly2004.html
double y = year+((month-1)*30.5+day/31*30.5)/366; double y = year+((month-1)*30.5+day/31.*30.5)/366;
// set the default value for Delta T // set the default value for Delta T
double u = (y-1820)/100.; double u = (y-1820)/100.;
double r = (-20 + 32 * u * u); double r = (-20 + 32 * u * u);
if (y < -500) if (y < -500)
{ {
// values are equal to defaults! // values are equal to defaults!
// u = (y-1820)/100.; // u = (y-1820)/100.;
// r = (-20 + 32 * u * u); // r = (-20 + 32 * u * u);
skipping to change at line 1347 skipping to change at line 1378
else if (y < 1961) else if (y < 1961)
{ {
double t = y - 1950; double t = y - 1950;
//r = (29.07 + 0.407*t - std::pow(t,2)/233 + std::pow(t,3) / 2547); //r = (29.07 + 0.407*t - std::pow(t,2)/233 + std::pow(t,3) / 2547);
r = ((t/2547.0 -1.0/233.0)*t + 0.407)*t +29.07; r = ((t/2547.0 -1.0/233.0)*t + 0.407)*t +29.07;
} }
else if (y < 1986) else if (y < 1986)
{ {
double t = y - 1975; double t = y - 1975;
//r = (45.45 + 1.067*t - std::pow(t,2)/260 - std::pow(t,3) / 718); //r = (45.45 + 1.067*t - std::pow(t,2)/260 - std::pow(t,3) / 718);
r = ((-t/718.0 +1/260.0)*t + 1.067)*t +45.45; r = ((-t/718.0 -1/260.0)*t + 1.067)*t +45.45;
} }
else if (y < 2005) else if (y < 2005)
{ {
double t = y - 2000; double t = y - 2000;
//r = (63.86 + 0.3345 * t - 0.060374 * std::pow(t,2) + 0.001 7275 * std::pow(t,3) + 0.000651814 * std::pow(t,4) + 0.00002373599 * std::p ow(t,5)); //r = (63.86 + 0.3345 * t - 0.060374 * std::pow(t,2) + 0.001 7275 * std::pow(t,3) + 0.000651814 * std::pow(t,4) + 0.00002373599 * std::p ow(t,5));
r = ((((0.00002373599*t + 0.000651814)*t + 0.0017275)*t - 0. 060374)*t + 0.3345)*t +63.86; r = ((((0.00002373599*t + 0.000651814)*t + 0.0017275)*t - 0. 060374)*t + 0.3345)*t +63.86;
} }
else if (y < 2050) else if (y < 2050)
{ {
double t = y - 2000; double t = y - 2000;
skipping to change at line 1375 skipping to change at line 1406
r -= 0.5628*(2150.0-y); r -= 0.5628*(2150.0-y);
} }
return r; return r;
} }
// Implementation of algorithm by Schoch (1931) for DeltaT computation // Implementation of algorithm by Schoch (1931) for DeltaT computation
double getDeltaTBySchoch(const double jDay) double getDeltaTBySchoch(const double jDay)
{ {
double u=(jDay-2378496.0)/36525.0; // (1800-jan-0.5) double u=(jDay-2378496.0)/36525.0; // (1800-jan-0.5)
return -36.28 + 36.28*std::pow(u,2); return -36.28 + 36.28*u*u;
} }
// Implementation of algorithm by Clemence (1948) for DeltaT computation // Implementation of algorithm by Clemence (1948) for DeltaT computation
double getDeltaTByClemence(const double jDay) double getDeltaTByClemence(const double jDay)
{ {
double u=(jDay-2415020.0)/36525.0; // (1900-jan-0.5) double u=(jDay-2415020.0)/36525.0; // (1900-jan-0.5)
return +8.72 + 26.75*u + 11.22*std::pow(u,2); return +8.72 + 26.75*u + 11.22*u*u;
} }
// Implementation of algorithm by IAU (1952) for DeltaT computation // Implementation of algorithm by IAU (1952) for DeltaT computation
double getDeltaTByIAU(const double jDay) double getDeltaTByIAU(const double jDay)
{ {
double u=(jDay-2415020.0)/36525.0; // (1900-jan-0.5) double u=(jDay-2415020.0)/36525.0; // (1900-jan-0.5)
// TODO: Calculate Moon's longitude fluctuation // TODO: Calculate Moon's longitude fluctuation
return (29.950*u +72.318)*u +24.349 /* + 1.82144*b */ ; return (29.950*u +72.318)*u +24.349 /* + 1.82144*b */ ;
} }
skipping to change at line 1443 skipping to change at line 1474
double u=(jDay-2415020.0)/36525.0; // (1900-jan-0.5) double u=(jDay-2415020.0)/36525.0; // (1900-jan-0.5)
double deltaT=(((((((((((-0.089491*u -0.117389)*u + 0.185489)*u + 0. 247433)*u - 0.159732)*u - 0.200097)*u + 0.075456)*u double deltaT=(((((((((((-0.089491*u -0.117389)*u + 0.185489)*u + 0. 247433)*u - 0.159732)*u - 0.200097)*u + 0.075456)*u
+ 0.076929)*u - 0.020446)*u - 0.013867)*u + 0.003081 )*u + 0.001233)*u -0.000029; + 0.076929)*u - 0.020446)*u - 0.013867)*u + 0.003081 )*u + 0.001233)*u -0.000029;
return deltaT * 86400.0; return deltaT * 86400.0;
} }
// Implementation of algorithm by Morrison & Stephenson (1982) for DeltaT c omputation // Implementation of algorithm by Morrison & Stephenson (1982) for DeltaT c omputation
double getDeltaTByMorrisonStephenson1982(const double jDay) double getDeltaTByMorrisonStephenson1982(const double jDay)
{ {
double u=(jDay-2382148.0)/36525.0; // (1810-jan-0.5) double u=(jDay-2382148.0)/36525.0; // (1810-jan-0.5)
return -15.0+32.50*std::pow(u,2); return -15.0+32.50*u*u;
} }
// Implementation of algorithm by Stephenson & Morrison (1984) for DeltaT c omputation // Implementation of algorithm by Stephenson & Morrison (1984) for DeltaT c omputation
double getDeltaTByStephensonMorrison1984(const double jDay) double getDeltaTByStephensonMorrison1984(const double jDay)
{ {
int year, month, day; int year, month, day;
double deltaT = 0.; double deltaT = 0.;
getDateFromJulianDay(jDay, &year, &month, &day); getDateFromJulianDay(jDay, &year, &month, &day);
double yeardec=year+((month-1)*30.5+day/31*30.5)/366; double u = (getDecYear(year, month, day)-1800)/100;
double u = (yeardec-1800)/100;
if (-391 < year && year <= 948) if (-391 < year && year <= 948)
deltaT = (44.3*u +320.0)*u +1360.0; deltaT = (44.3*u +320.0)*u +1360.0;
if (948 < year && year <= 1600) if (948 < year && year <= 1600)
deltaT = 25.5*u*u; deltaT = 25.5*u*u;
return deltaT; return deltaT;
} }
// Implementation of algorithm by Stephenson & Morrison (1995) for DeltaT c omputation // Implementation of algorithm by Stephenson & Morrison (1995) for DeltaT c omputation
skipping to change at line 1479 skipping to change at line 1509
} }
// Implementation of algorithm by Stephenson & Houlden (1986) for DeltaT co mputation // Implementation of algorithm by Stephenson & Houlden (1986) for DeltaT co mputation
double getDeltaTByStephensonHoulden(const double jDay) double getDeltaTByStephensonHoulden(const double jDay)
{ {
int year, month, day; int year, month, day;
double u; double u;
double deltaT = 0.; double deltaT = 0.;
getDateFromJulianDay(jDay, &year, &month, &day); getDateFromJulianDay(jDay, &year, &month, &day);
double yeardec=year+((month-1)*30.5+day/31*30.5)/366; double yeardec=getDecYear(year, month, day);
if (year <= 948) if (year <= 948)
{ {
u = (yeardec-948)/100; u = (yeardec-948)/100;
deltaT = (46.5*u -405.0)*u + 1830.0; deltaT = (46.5*u -405.0)*u + 1830.0;
} }
if (948 < year && year <= 1600) if (948 < year && year <= 1600)
{ {
u = (yeardec-1850)/100; u = (yeardec-1850)/100;
deltaT = 22.5*u*u; deltaT = 22.5*u*u;
skipping to change at line 1569 skipping to change at line 1599
{ {
double u=(jDay-2385800.0)/36525.0; // (1820-jan-0.5) double u=(jDay-2385800.0)/36525.0; // (1820-jan-0.5)
return -20.0 + 32.0 * u*u; return -20.0 + 32.0 * u*u;
} }
// Implementation of algorithm by Reijs (2006) for DeltaT computation // Implementation of algorithm by Reijs (2006) for DeltaT computation
double getDeltaTByReijs(const double jDay) double getDeltaTByReijs(const double jDay)
{ {
double OffSetYear = (2385800.0 - jDay)/365.25; double OffSetYear = (2385800.0 - jDay)/365.25;
return ((1.8 * std::pow(OffSetYear,2)/200 + 1443*3.76/(2*M_PI)*(std: :cos(2*M_PI*OffSetYear/1443)-1))*365.25)/1000; return ((1.8 * OffSetYear*OffSetYear/200 + 1443*3.76/(2*M_PI)*(std:: cos(2*M_PI*OffSetYear/1443)-1))*365.25)/1000;
} }
// Implementation of algorithm by Chapront, Chapront-Touze & Francou (1997) & Meeus (1998) for DeltaT computation // Implementation of algorithm by Chapront, Chapront-Touze & Francou (1997) & Meeus (1998) for DeltaT computation
// 191 values in tenths of second for interpolation table 1620..2000, every 2nd year. // 191 values in tenths of second for interpolation table 1620..2000, every 2nd year.
static const int MeeusDeltaTTable[] = { 1210, 1120, 1030, 950, 880, 820, 77 0, 720, 680, 630, 600, 560, 530, 510, 480, static const int MeeusDeltaTTable[] = { 1210, 1120, 1030, 950, 880, 820, 77 0, 720, 680, 630, 600, 560, 530, 510, 480,
460, 440, 420, 400, 380, 350, 330, 3 10, 290, 260, 240, 220, 200, 180, 160, 140, 120, 110, 100, 90, 80, 70, 70, 70, 70, // before 1700 460, 440, 420, 400, 380, 350, 330, 3 10, 290, 260, 240, 220, 200, 180, 160, 140, 120, 110, 100, 90, 80, 70, 70, 70, 70, // before 1700
70, 70, 80, 80, 90, 90, 90, 9 0, 90, 100, 100, 100, 100, 100, 100, 100, 100, 110, 110, 110, 110, 110, 12 0, 120, 120, // before 1750 70, 70, 80, 80, 90, 90, 90, 9 0, 90, 100, 100, 100, 100, 100, 100, 100, 100, 110, 110, 110, 110, 110, 12 0, 120, 120, // before 1750
120, 130, 130, 130, 140, 140, 140, 1 40, 150, 150, 150, 150, 150, 160, 160, 160, 160, 160, 160, 160, 160, 150, 1 50, 140, 130, // before 1800 120, 130, 130, 130, 140, 140, 140, 1 40, 150, 150, 150, 150, 150, 160, 160, 160, 160, 160, 160, 160, 160, 150, 1 50, 140, 130, // before 1800
131, 125, 122, 120, 120, 120, 120, 1 20, 120, 119, 116, 110, 102, 92, 82, 71, 62, 56, 54, 53, 54, 56, 59, 62, 65, // before 1850 131, 125, 122, 120, 120, 120, 120, 1 20, 120, 119, 116, 110, 102, 92, 82, 71, 62, 56, 54, 53, 54, 56, 59, 62, 65, // before 1850
68, 71, 73, 75, 76, 77, 73, 6 2, 52, 27, 14, -12, -28, -38, -48, -55, -53, -56, -57, -59, -60, -63, -6 5, -62, -47, // before 1900 68, 71, 73, 75, 76, 77, 73, 6 2, 52, 27, 14, -12, -28, -38, -48, -55, -53, -56, -57, -59, -60, -63, -6 5, -62, -47, // before 1900
skipping to change at line 1603 skipping to change at line 1633
else if (year<1620) else if (year<1620)
deltaT=(25.3*u20 + 102.0)*u20+102.0; deltaT=(25.3*u20 + 102.0)*u20+102.0;
// The next terms are the short approximative formulae. But Meeus gi ves exact values, table, etc. for 1620..2000. // The next terms are the short approximative formulae. But Meeus gi ves exact values, table, etc. for 1620..2000.
//else if (year<1900) //else if (year<1900)
// deltaT= ((((((((( 123563.95*u19 + 727058.63)*u19 + 1818961 .41)*u19 + 2513807.78)*u19 + 2087298.89)*u19 + 1061660.75)*u19 + // deltaT= ((((((((( 123563.95*u19 + 727058.63)*u19 + 1818961 .41)*u19 + 2513807.78)*u19 + 2087298.89)*u19 + 1061660.75)*u19 +
// 324011.78)*u19 + 56282.84)*u19 + 5218.61)*u19 + 228.95)*u19 - 2.50; // 324011.78)*u19 + 56282.84)*u19 + 5218.61)*u19 + 228.95)*u19 - 2.50;
//else if (year<1997) //else if (year<1997)
// deltaT= (((((((( 58353.42*u19 -232424.66)*u19 +372919.88)* u19 - 303191.19)*u19 + 124906.15)*u19 - 18756.33)*u19 - 2637.80)*u19 + 815. 20)*u19 + 87.24)*u19 - 2.44; // deltaT= (((((((( 58353.42*u19 -232424.66)*u19 +372919.88)* u19 - 303191.19)*u19 + 124906.15)*u19 - 18756.33)*u19 - 2637.80)*u19 + 815. 20)*u19 + 87.24)*u19 - 2.44;
else if (year <2000) else if (year <2000)
{ {
double yeardec=year+((month-1)*30.5+day/31*30.5)/366; double yeardec=getDecYear(year, month, day);
int pos=(year-1620)/2; // this is a deliberate integer divis ion! 2->1, 3->1, 4->2, 5->2 etc. int pos=(year-1620)/2; // this is a deliberate integer divis ion! 2->1, 3->1, 4->2, 5->2 etc.
deltaT= MeeusDeltaTTable[pos]+ (yeardec-(2*pos+1620))*0.5 * (MeeusDeltaTTable[pos+1]-MeeusDeltaTTable[pos]); deltaT= MeeusDeltaTTable[pos]+ (yeardec-(2*pos+1620))*0.5 * (MeeusDeltaTTable[pos+1]-MeeusDeltaTTable[pos]);
deltaT /= 10.0; deltaT /= 10.0;
} }
else if (year<2100) else if (year<2100)
deltaT=(25.3*u20 + 102.0)*u20+102.0 + 0.37*(year-2100); deltaT=(25.3*u20 + 102.0)*u20+102.0 + 0.37*(year-2100);
else // year > 2100 else // year > 2100
deltaT=(25.3*u20 + 102.0)*u20+102.0; deltaT=(25.3*u20 + 102.0)*u20+102.0;
return deltaT; return deltaT;
} }
skipping to change at line 1842 skipping to change at line 1872
else if (year <= 2007) else if (year <= 2007)
{ {
// revised 2013 per email // revised 2013 per email
u = 0.265 + ub; u = 0.265 + ub;
deltaT = (((1314.759 * u - 296.018) * u - 101.898) * u + 88. 659) * u + 49.997; deltaT = (((1314.759 * u - 296.018) * u - 101.898) * u + 88. 659) * u + 49.997;
} }
return deltaT; return deltaT;
} }
// Implementation of polinomial approximation of time period 1620-2013 for
DeltaT by M. Khalid, Mariam Sultana and Faheem Zaidi (2014).
double getDeltaTByKhalidSultanaZaidi(const double jDay)
{
int year, month, day;
getDateFromJulianDay(jDay, &year, &month, &day);
double k, a0, a1, a2, a3, a4;
if (year>=1620 && year<=1672)
{
k = 3.670; a0 = 76.541; a1 = -253.532; a2 = 695.901; a3 = -1
256.982; a4 = 627.152;
}
else if (year>=1673 && year<=1729)
{
k = 3.120; a0 = 10.872; a1 = -40.744; a2 = 236.890; a3 = -35
1.537; a4 = 36.612;
}
else if (year>=1730 && year<=1797)
{
k = 2.495; a0 = 13.480; a1 = 13.075; a2 = 8.635; a3 = -3.307
; a4 = -128.294;
}
else if (year>=1798 && year<=1843)
{
k = 1.925; a0 = 12.584; a1 = 1.929; a2 = 60.896; a3 = -1432.
216; a4 = 3129.071;
}
else if (year>=1844 && year<=1877)
{
k = 1.525; a0 = 6.364; a1 = 11.004; a2 = 407.776; a3 = -4168
.394; a4 = 7561.686;
}
else if (year>=1878 && year<=1904)
{
k = 1.220; a0 = -5.058; a1 = -1.701; a2 = -46.403; a3 = -866
.171; a4 = 5917.585;
}
else if (year>=1905 && year<=1945)
{
k = 0.880; a0 = 13.392; a1 = 128.592; a2 = -279.165; a3 = -1
282.050; a4 = 4039.490;
}
else if (year>=1946 && year<=1989)
{
k = 0.455; a0 = 30.782; a1 = 34.348; a2 = 46.452; a3 = 1295.
550; a4 = -3210.913;
}
else if (year>=1990 && year<=2013)
{
k = 0.115; a0 = 55.281; a1 = 91.248; a2 = 87.202; a3 = -3092
.565; a4 = 8255.422;
}
else
{
k = 0.0; a0 = 0.0; a1 = 0.0; a2 = 0.0; a3 = 0.0; a4 = 0.0;
}
double u = k + (year - 2000)/100;
return (((a4*u + a3)*u + a2)*u + a1)*u + a0;
}
double getMoonSecularAcceleration(const double jDay, const double nd) double getMoonSecularAcceleration(const double jDay, const double nd)
{ {
int year, month, day; int year, month, day;
getDateFromJulianDay(jDay, &year, &month, &day); getDateFromJulianDay(jDay, &year, &month, &day);
double yeardec=year+((month-1)*30.5+day/31*30.5)/366.0; double t = (getDecYear(year, month, day)-1955.5)/100.0;
double t = (yeardec-1955.5)/100.0;
// n.dot for secular acceleration of the Moon in ELP2000-82B // n.dot for secular acceleration of the Moon in ELP2000-82B
// have value -23.8946 "/cy/cy // have value -23.8946 "/cy/cy
return -0.91072 * (-23.8946 + qAbs(nd))*t*t; return -0.91072 * (-23.8946 + qAbs(nd))*t*t;
} }
double getDeltaTStandardError(const double jDay) double getDeltaTStandardError(const double jDay)
{ {
int year, month, day; int year, month, day;
getDateFromJulianDay(jDay, &year, &month, &day); getDateFromJulianDay(jDay, &year, &month, &day);
//double yeardec=year+((month-1)*30.5+day/31*30.5)/366;
double sigma = -1.; double sigma = -1.;
if (-1000 <= year && year <= 1600) if (-1000 <= year && year <= 1600)
{ {
double cDiff1820= (jDay-2385800.0)/36525.0; // 1820.0=182 0-jan-0.5=2385800.0 double cDiff1820= (jDay-2385800.0)/36525.0; // 1820.0=182 0-jan-0.5=2385800.0
// sigma = std::pow((yeardec-1820.0)/100,2); // sigma(DeltaT ) = 0.8*u^2 // sigma = std::pow((yeardec-1820.0)/100,2); // sigma(DeltaT ) = 0.8*u^2
sigma = 0.8 * cDiff1820 * cDiff1820; sigma = 0.8 * cDiff1820 * cDiff1820;
} }
return sigma; return sigma;
} }
skipping to change at line 1969 skipping to change at line 2049
*cos_sin++ = sin(minAngle); *cos_sin++ = sin(minAngle);
for (int i=0; i<segments; ++i) // we cannot mirror this, it may be u nequal. for (int i=0; i<segments; ++i) // we cannot mirror this, it may be u nequal.
{ // efficient computation, avoid expensive trig functions by use of the addition theorem. { // efficient computation, avoid expensive trig functions by use of the addition theorem.
cos_sin[0] = cos_sin[-2]*c - cos_sin[-1]*s; cos_sin[0] = cos_sin[-2]*c - cos_sin[-1]*s;
cos_sin[1] = cos_sin[-2]*s + cos_sin[-1]*c; cos_sin[1] = cos_sin[-2]*s + cos_sin[-1]*c;
cos_sin += 2; cos_sin += 2;
} }
return cos_sin_rho; return cos_sin_rho;
} }
double getDecYear(const int year, const int month, const int day)
{
return year+((month-1)*30.5+day/31.*30.5)/366;
}
//! Uncompress gzip or zlib compressed data. //! Uncompress gzip or zlib compressed data.
QByteArray uncompress(const QByteArray& data) QByteArray uncompress(const QByteArray& data)
{ {
if (data.size() <= 4) if (data.size() <= 4)
return QByteArray(); return QByteArray();
static const int CHUNK = 1024;
QByteArray buffer(CHUNK, 0); //needed for const-correctness, no deep copy performed
QByteArray dataNonConst(data);
QBuffer buf(&dataNonConst);
buf.open(QIODevice::ReadOnly);
return uncompress(buf);
}
//! Uncompress (gzip/zlib) data from this QIODevice, which must be open and
readable.
//! @param device The device to read from, must already be opened with an O
penMode supporting reading
//! @param maxBytes The max. amount of bytes to read from the device, or -1
to read until EOF. Note that it
//! always stops when inflate() returns Z_STREAM_END. Positive values can b
e used for interleaving compressed data
//! with other data.
QByteArray uncompress(QIODevice& device, qint64 maxBytes)
{
// this is a basic zlib decompression routine, similar to:
// http://zlib.net/zlib_how.html
// buffer size 256k, zlib recommended size
static const int CHUNK = 262144;
QByteArray readBuffer(CHUNK, 0);
QByteArray inflateBuffer(CHUNK, 0);
QByteArray out; QByteArray out;
// zlib stream
z_stream strm; z_stream strm;
strm.zalloc = Z_NULL; strm.zalloc = Z_NULL;
strm.zfree = Z_NULL; strm.zfree = Z_NULL;
strm.opaque = Z_NULL; strm.opaque = Z_NULL;
strm.avail_in = data.size(); strm.avail_in = Z_NULL;
strm.next_in = (Bytef*)(data.data()); strm.next_in = Z_NULL;
// the amount of bytes already read from the device
qint64 bytesRead = 0;
// initialize zlib
// 15 + 32 for gzip automatic header detection. // 15 + 32 for gzip automatic header detection.
int ret = inflateInit2(&strm, 15 + 32); int ret = inflateInit2(&strm, 15 + 32);
if (ret != Z_OK) return QByteArray(); if (ret != Z_OK)
{
qWarning()<<"zlib init error ("<<ret<<"), can't uncompress";
if(strm.msg)
qWarning()<<"zlib message: "<<QString(strm.msg);
return QByteArray();
}
//zlib double loop - one for reading from file, one for inflating
do do
{ {
strm.avail_out = CHUNK; qint64 bytesToRead = CHUNK;
strm.next_out = (Bytef*)(buffer.data()); if(maxBytes>=0)
ret = inflate(&strm, Z_NO_FLUSH);
Q_ASSERT(ret != Z_STREAM_ERROR);
if (ret < 0)
{ {
out.clear(); //check if we reach the desired limit with the next
read
bytesToRead = qMin((qint64)CHUNK,maxBytes-bytesRead)
;
}
if(bytesToRead==0)
break; break;
//perform read from device
qint64 read = device.read(readBuffer.data(), bytesToRead);
if (read<0)
{
qWarning()<<"Error while reading from device";
inflateEnd(&strm);
return QByteArray();
} }
out.append(buffer.data(), CHUNK - strm.avail_out);
bytesRead += read;
strm.next_in = reinterpret_cast<Bytef*>(readBuffer.data());
strm.avail_in = read;
if(read==0)
break;
//inflate loop
do
{
strm.avail_out = CHUNK;
strm.next_out = reinterpret_cast<Bytef*>(inflateBuff
er.data());
ret = inflate(&strm,Z_NO_FLUSH);
Q_ASSERT(ret != Z_STREAM_ERROR); // must never happe
n, indicates a programming error
if(ret < 0 || ret == Z_NEED_DICT)
{
qWarning()<<"zlib inflate error ("<<ret<<"),
can't uncompress";
if(strm.msg)
qWarning()<<"zlib message: "<<QStrin
g(strm.msg);
inflateEnd(&strm);
return QByteArray();
}
out.append(inflateBuffer.constData(), CHUNK - strm.a
vail_out);
}while(strm.avail_out == 0); //if zlib has more data for us,
repeat
}while(ret!=Z_STREAM_END);
// close zlib
inflateEnd(&strm);
if(ret!=Z_STREAM_END)
{
qWarning()<<"Premature end of compressed stream";
if(strm.msg)
qWarning()<<"zlib message: "<<QString(strm.msg);
return QByteArray();
} }
while (strm.avail_out == 0);
inflateEnd(&strm);
return out; return out;
} }
} // end of the StelUtils namespace } // end of the StelUtils namespace
 End of changes. 29 change blocks. 
29 lines changed or deleted 217 lines changed or added

This html diff was produced by rfcdiff 1.41. The latest version is available from http://tools.ietf.org/tools/rfcdiff/