PLplot
5.10.0
|
00001 // 00002 // This software originally contributed under the LGPL in January 2009 to 00003 // PLplot by the 00004 // Cluster Science Centre 00005 // QSAS team, 00006 // Imperial College, London 00007 // Copyright (C) 2009 Imperial College, London 00008 // 00009 // This file is part of PLplot. 00010 // 00011 // PLplot is free software; you can redistribute it and/or modify 00012 // it under the terms of the GNU Library General Public License as published 00013 // by the Free Software Foundation; either version 2 of the License, or 00014 // (at your option) any later version. 00015 // 00016 // PLplot is distributed in the hope that it will be useful, 00017 // but WITHOUT ANY WARRANTY; without even the implied warranty of 00018 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 00019 // GNU Library General Public License for more details. 00020 // 00021 // You should have received a copy of the GNU Library General Public License 00022 // along with PLplot; if not, write to the Free Software 00023 // Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA 00024 // 00025 // 00026 00027 // MJD measures from the start of 17 Nov 1858 00028 00029 // These utilities use the Gregorian calendar after 4 Oct 1582 (Julian) i.e. from 15 Oct 1582 Gregoriam 00030 // Note C libraries use Gregorian only from 14 Sept 1752 00031 // More detailed discussion can be found at http://aa.usno.navy.mil/data/docs/JulianDate.php 00032 // These routines have been compared with the results of the US Naval Observatory online converter. 00033 // Modified Julian Date (MJD) = Julian Date (JD) - 2400000.5 00034 // 00035 // In all routines, specifying a day, hour, minute or second field greater than would be valid is 00036 // handled with modulo arithmetic and safe. 00037 // Thus 2006-12-32 00:62:00.0 will safely, and correctly, be treated as 2007-01-01 01:02:00.0 00038 // 00039 // 00040 00041 #include "qsastime_extra.h" 00042 00043 static double MJDtoJD = 2400000.5; 00044 static double SecInDay = 86400; // we ignore leap seconds 00045 00046 int setFromISOstring( const char* ISOstring, MJDtime *MJD, int forceJulian ) 00047 { 00048 double seconds; 00049 int y, m, d, h, min; 00050 int startAt = 0; 00051 int len = strlen( ISOstring ); 00052 00053 // ISO is "1995-01-23 02:33:17.235" or "1995-01-23T02:33:17.235Z" 00054 00055 // parse off year 00056 00057 y = strtol( &( ISOstring[startAt] ), NULL, 10 ); 00058 if ( ISOstring[startAt] == '-' || ISOstring[startAt] == '+' ) 00059 startAt++; 00060 startAt += 5; 00061 if ( startAt > len ) 00062 return 1; 00063 00064 m = strtol( &( ISOstring[startAt] ), NULL, 10 ); 00065 startAt += 3; 00066 if ( startAt > len ) 00067 return 1; 00068 00069 d = strtol( &( ISOstring[startAt] ), NULL, 10 ); 00070 startAt += 3; 00071 if ( startAt > len ) 00072 return 1; 00073 00074 h = strtol( &( ISOstring[startAt] ), NULL, 10 ); 00075 startAt += 3; 00076 if ( startAt > len ) 00077 return 1; 00078 00079 min = strtol( &( ISOstring[startAt] ), NULL, 10 ); 00080 startAt += 3; 00081 if ( startAt > len ) 00082 return 1; 00083 00084 seconds = strtod( &( ISOstring[startAt] ), NULL ); 00085 setFromUT( y, m - 1, d, h, min, seconds, MJD, forceJulian ); 00086 00087 return 0; 00088 } 00089 00090 00091 void setFromDOY( int year, int doy, int hour, int min, double sec, MJDtime *MJD, int forceJulian ) 00092 { 00093 // Set from Day Of Year format 00094 00095 // convert Gregorian date plus time to MJD 00096 // MJD measures from the start of 17 Nov 1858 00097 00098 // the int flag forceJulian forces use of Julian calendar whatever the year 00099 // default is to use Gregorian after 4 Oct 1582 (Julian) i.e. from 15 Oct 1582 Gregorian 00100 // Note C libraries use Gregorian only from 14 Sept 1752 onwards 00101 00102 int leaps, lastyear, extraDays; 00103 00104 // N.B. There were known problems (both for the Julian and Gregorian 00105 // cases) with the following leap year logic that were completely fixed 00106 // in qsastime.c, but I (AWI) am not going to bother with these fixups 00107 // here since this code only used for a specific test routine for limited 00108 // date range and not for anything general. 00109 if ( forceJulian && year <= 0 ) 00110 { 00111 // count leap years on Julian Calendar 00112 // MJD for Jan 1 0000 (correctly Jan 01, BCE 1) is - 678943, count from there 00113 // negative CE (AD) years convert to BCE (BC) as BCE = 1 - CE, e.g. 2 BCE = -1 CE 00114 00115 leaps = ( year - 4 ) / 4; // (note leaps is negative here and year 0 (1 BCE) was a leap year 00116 MJD->base_day = year * 365 + leaps + doy - 678943; 00117 } 00118 else if ( forceJulian ) 00119 { 00120 // count leap years on Julian Calendar 00121 // MJD for Jan 1 0000 (correctly Jan 01, BCE 1) is - 678943, count from there 00122 00123 leaps = ( year - 1 ) / 4; 00124 MJD->base_day = year * 365 + leaps + doy - 678943; 00125 } 00126 else 00127 { 00128 // count leap years Gregorian Calendar - modern dates 00129 // Algorithm below for 17 Nov 1858 (0 MJD) gives 00130 // leaps = 450 and hence base_day of 678941, so subtract it to give MJD day 00131 00132 lastyear = year - 1; 00133 leaps = lastyear / 4 - lastyear / 100 + lastyear / 400; 00134 MJD->base_day = year * 365 + leaps + doy - 678941; 00135 } 00136 00137 MJD->time_sec = sec + ( (double) min + (double) hour * 60. ) * 60.; 00138 00139 if ( MJD->time_sec >= SecInDay ) 00140 { 00141 extraDays = (int) ( MJD->time_sec / SecInDay ); 00142 MJD->base_day += extraDays; 00143 MJD->time_sec -= extraDays * SecInDay; 00144 } 00145 00146 return; 00147 } 00148 00149 00150 void setFromBCE( int yearBCE, int month, int day, int hour, int min, double sec, MJDtime *MJD, int forceJulian ) 00151 { 00152 // utility to allow user to input dates BCE (BC) 00153 00154 int year = 1 - yearBCE; 00155 setFromUT( year, month, day, hour, min, sec, MJD, forceJulian ); 00156 } 00157 00158 void setFromMJD( double ModifiedJulianDate, MJDtime *MJD ) 00159 { 00160 // convert MJD double into MJD structure 00161 MJD->base_day = (int) ModifiedJulianDate; 00162 MJD->time_sec = ( ModifiedJulianDate - MJD->base_day ) * SecInDay; 00163 } 00164 00165 void setFromJD( double JulianDate, MJDtime *MJD ) 00166 { 00167 // break JD double into MJD based structure 00168 // Note Julian Day starts Noon, so convert to MJD first 00169 00170 MJD->base_day = (int) ( JulianDate - MJDtoJD ); 00171 MJD->time_sec = ( JulianDate - MJDtoJD - (double) MJD->base_day ) * SecInDay; 00172 } 00173 00174 void setFromCDFepoch( double cdfepoch, MJDtime *MJD ) 00175 { 00176 // convert cdf epoch double into MJD structure 00177 // Note that cdfepoch is msec from 0 AD on the Gregorian calendar 00178 00179 double seconds = cdfepoch * 0.001; 00180 00181 MJD->base_day = (int) ( seconds / 86400.0 ); 00182 MJD->time_sec = seconds - MJD->base_day * SecInDay; 00183 MJD->base_day -= 678941; 00184 } 00185 00186 double getCDFepoch( MJDtime *MJD ) 00187 { 00188 // convert MJD structure into cdf epoch double 00189 // Note that cdfepoch is msec from 0 AD on the Gregorian Calendar 00190 00191 int days = MJD->base_day + 678941; 00192 double seconds = days * SecInDay + MJD->time_sec; 00193 return seconds * 1000.; 00194 } 00195 00196 double getMJD( MJDtime *MJD ) 00197 { 00198 // Return MJD as a double 00199 return (double) MJD->base_day + MJD->time_sec / SecInDay; 00200 } 00201 00202 double getJD( MJDtime *MJD ) 00203 { 00204 // Return JD as a double 00205 double JD = getMJD( MJD ) + MJDtoJD; 00206 return JD; 00207 } 00208 00209 double getDiffDays( MJDtime *MJD1, MJDtime *MJD2 ) 00210 { 00211 // Return difference MJD1 - MJD2 in days as a double 00212 double diff = (double) ( MJD1->base_day - MJD2->base_day ) + ( MJD1->time_sec - MJD2->time_sec ) / SecInDay; 00213 return diff; 00214 } 00215 00216 double getDiffSecs( MJDtime *MJD1, MJDtime *MJD2 ) 00217 { 00218 // Return difference MJD1 - MJD2 in seconds as a double 00219 double diff = (double) ( MJD1->base_day - MJD2->base_day ) * SecInDay + ( MJD1->time_sec - MJD2->time_sec ); 00220 return diff; 00221 } 00222 00223 const char * getISOString( MJDtime* MJD, int delim, int forceJulian ) 00224 { 00225 // ISO time string for UTC 00226 // uses default behaviour for Julian/Gregorian switch over 00227 //** 00228 // Warning getISOString is not thread safe 00229 // as it writes to a static variable DateTime 00230 //* 00231 00232 static char DateTime[50]; 00233 int y, m, d, hour, min; 00234 int sec1, ysign; 00235 double sec; 00236 int slen; 00237 char * ptr; 00238 00239 breakDownMJD( &y, &m, &d, &hour, &min, &sec, MJD, forceJulian ); 00240 00241 if ( y < 0 ) 00242 { 00243 ysign = 1; 00244 y = -y; 00245 } 00246 else 00247 ysign = 0; 00248 00249 sec1 = (int) sec / 10; 00250 sec -= (double) sec1 * 10; 00251 00252 if ( delim == 1 ) 00253 { 00254 if ( ysign == 0 ) 00255 sprintf( DateTime, "%04d-%02d-%02dT%02d:%02d:%01d%-11.10f", y, m + 1, d, hour, min, sec1, sec ); 00256 else 00257 sprintf( DateTime, "-%04d-%02d-%02dT%02d:%02d:%01d%-11.10f", y, m + 1, d, hour, min, sec1, sec ); 00258 00259 // remove trailing white space 00260 while ( ( ptr = strrchr( &( DateTime[0] ), ' ' ) ) != NULL ) 00261 ptr[0] = '\0'; 00262 strcat( &( DateTime[0] ), "Z" ); 00263 } 00264 else 00265 { 00266 if ( ysign == 0 ) 00267 sprintf( DateTime, "%04d-%02d-%02d %02d:%02d:%01d%-11.10f", y, m + 1, d, hour, min, sec1, sec ); 00268 else 00269 sprintf( DateTime, "-%04d-%02d-%02d %02d:%02d:%01d%-11.10f", y, m + 1, d, hour, min, sec1, sec ); 00270 00271 // remove trailing white space 00272 slen = strlen( DateTime ) - 1; 00273 while ( DateTime[slen] == ' ' ) 00274 { 00275 DateTime[slen] = '\0'; 00276 slen--; 00277 } 00278 } 00279 return &( DateTime[0] ); 00280 }