external/io_grib1/WGRIB/PDS_date.c

   1 #include <stdio.h>
   2 #include <stdlib.h>
   3 #include <stddef.h>
   4 #include <string.h>
   5 #include "pds4.h"
   6 #include "grib.h"
   7
   8 /*
   9  * PDS_date.c  v1.2 wesley ebisuzaki
  10  *
  11  * prints a string with a date code
  12  *
  13  * PDS_date(pds,option, v_time)
  14  *   options=0  .. 2 digit year
  15  *   options=1  .. 4 digit year
  16  *
  17  *   v_time=0   .. initial time
  18  *   v_time=1   .. verification time
  19  *
  20  * assumption: P1 and P2 are unsigned integers (not clear from doc)
  21  *
  22  * v1.2 years that are multiple of 400 are leap years, not 500
  23  * v1.2.1  make the change to the source code for v1.2
  24  * v1.2.2  add 3/6/12 hour forecast time units
  25  */
  26
  27 static int msg_count = 0;
  28 extern int minute;
  29
  30 int PDS_date(unsigned char *pds, int option, int v_time) {
  31
  32     int year, month, day, hour, min;
  33
  34     if (v_time == 0) {
  35         year = PDS_Year4(pds);
  36         month = PDS_Month(pds);
  37         day  = PDS_Day(pds);
  38         hour = PDS_Hour(pds);
  39     }
  40     else {
  41         if (verf_time(pds, &year, &month, &day, &hour) != 0) {
  42             if (msg_count++ < 5) fprintf(stderr, "PDS_date: problem\n");
  43         }
  44     }
  45     min =  PDS_Minute(pds);
  46
  47     switch(option) {
  48         case 0:
  49             printf("%2.2d%2.2d%2.2d%2.2d", year % 100, month, day, hour);
  50             if (minute) printf("-%2.2d", min);
  51             break;
  52         case 1:
  53             printf("%4.4d%2.2d%2.2d%2.2d", year, month, day, hour);
  54             if (minute) printf("-%2.2d", min);
  55             break;
  56         default:
  57             fprintf(stderr,"missing code\n");
  58             exit(8);
  59     }
  60     return 0;
  61 }
  62
  63 #define  FEB29   (31+29)
  64 static int monthjday[12] = {
  65         0,31,59,90,120,151,181,212,243,273,304,334};
  66
  67 static int leap(int year) {
  68         if (year % 4 != 0) return 0;
  69         if (year % 100 != 0) return 1;
  70         return (year % 400 == 0);
  71 }
  72
  73
  74 int add_time(int *year, int *month, int *day, int *hour, int dtime, int unit) {
  75     int y, m, d, h, jday, i;
  76
  77     y = *year;
  78     m = *month;
  79     d = *day;
  80     h = *hour;
  81
  82     if (unit == YEAR) {
  83         *year = y + dtime;
  84         return 0;
  85     }
  86     if (unit == DECADE) {
  87         *year =  y + (10 * dtime);
  88         return 0;
  89     }
  90     if (unit == CENTURY) {
  91         *year = y + (100 * dtime);
  92         return 0;
  93     }
  94     if (unit == NORMAL) {
  95         *year = y + (30 * dtime);
  96         return 0;
  97     }
  98     if (unit == MONTH) {
  99         dtime += (m - 1);
 100         *year = y + (dtime / 12);
 101         *month = 1 + (dtime % 12);
 102         return 0;
 103     }
 104
 105     if (unit == SECOND) {
 106         dtime /= 60;
 107         unit = MINUTE;
 108     }
 109     if (unit == MINUTE) {
 110         dtime /= 60;
 111         unit = HOUR;
 112     }
 113
 114     if (unit == HOURS3) {
 115         dtime *= 3;
 116         unit = HOUR;
 117     }
 118     else if (unit == HOURS6) {
 119         dtime *= 6;
 120         unit = HOUR;
 121     }
 122     else if (unit == HOURS12) {
 123         dtime *= 12;
 124         unit = HOUR;
 125     }
 126
 127     if (unit == HOUR) {
 128         dtime += h;
 129         *hour = dtime % 24;
 130         dtime = dtime / 24;
 131         unit = DAY;
 132     }
 133
 134     /* this is the hard part */
 135
 136     if (unit == DAY) {
 137         /* set m and day to Jan 0, and readjust dtime */
 138         jday = d + monthjday[m-1];
 139         if (leap(y) && m > 2) jday++;
 140         dtime += jday;
 141
 142         /* 4 year chuncks */
 143         i = dtime / (4 * 365 + 1);
 144         if (i) {
 145             /* assume century years are leap */
 146             y = y + i*4;
 147             dtime -= i*(4 * 365 + 1);
 148             /* see if we have gone past feb 28, 1900, 2000, etc */
 149             if ((y - 1) / 100 != (*year-1) / 100) {
 150                 /* crossed the feb 28, xx00 */
 151                 /* correct for only one century mark */
 152                 if ((y / 100) % 4 != 0) dtime++;
 153             }
 154         }
 155
 156         /* one year chunks */
 157         while (dtime > 365 + leap(y)) {
 158             dtime -= (365 + leap(y));
 159             y++;
 160         }
 161
 162         /* calculate the month and day */
 163
 164         if (leap(y) && dtime == FEB29) {
 165             m = 2;
 166             d = 29;
 167         }
 168         else {
 169             if (leap(y) && dtime > FEB29) dtime--;
 170             for (i = 11; monthjday[i] >= dtime; --i);
 171             m = i + 1;
 172             d = dtime - monthjday[i];
 173         }
 174         *year = y;
 175         *month = m;
 176         *day = d;
 177         return 0;
 178    }
 179    fprintf(stderr,"add_time: undefined time unit %d\n", unit);
 180    return 1;
 181 }
 182
 183
 184 /*
 185  * verf_time:
 186  *
 187  * this routine returns the "verification" time
 188  * should have behavior similar to gribmap
 189  *
 190  */
 191
 192 int verf_time(unsigned char *pds, int *year, int *month, int *day, int *hour) {
 193     int tr, dtime, unit;
 194
 195     *year = PDS_Year4(pds);
 196     *month = PDS_Month(pds);
 197     *day  = PDS_Day(pds);
 198     *hour = PDS_Hour(pds);
 199
 200     /* find time increment */
 201
 202     dtime = PDS_P1(pds);
 203     tr = PDS_TimeRange(pds);
 204     unit = PDS_ForecastTimeUnit(pds);
 205
 206     if (tr == 10) dtime = PDS_P1(pds) * 256 + PDS_P2(pds);
 207     if (tr > 1 && tr < 6 ) dtime = PDS_P2(pds);
 208
 209     if (dtime == 0) return 0;
 210
 211     return add_time(year, month, day, hour, dtime, unit);
 212 }
 213