ViewVC Help
View File | Revision Log | Show Annotations | Download File | Root Listing
root/radiance/ray/src/gen/gensky.c
(Generate patch)

Comparing ray/src/gen/gensky.c (file contents):
Revision 2.22 by schorsch, Sun Nov 16 10:29:38 2003 UTC vs.
Revision 2.29 by greg, Thu Mar 30 20:19:05 2023 UTC

# Line 10 | Line 10 | static const char      RCSid[] = "$Id$";
10   *     3/26/86
11   */
12  
13 < #include  <stdio.h>
13 > #include  "rtio.h"
14   #include  <stdlib.h>
15 #include  <string.h>
15   #include  <math.h>
16   #include  <ctype.h>
17 <
17 > #include  "sun.h"
18   #include  "color.h"
19  
21 extern int jdate(int month, int day);
22 extern double stadj(int  jd);
23 extern double sdec(int  jd);
24 extern double salt(double sd, double st);
25 extern double sazi(double sd, double st);
26
20   #ifndef  PI
21   #define  PI             3.14159265358979323846
22   #endif
# Line 38 | Line 31 | extern double sazi(double sd, double st);
31   #define  overcast       ((skytype==S_OVER)|(skytype==S_UNIF))
32  
33   double  normsc();
41                                        /* sun calculation constants */
42 extern double  s_latitude;
43 extern double  s_longitude;
44 extern double  s_meridian;
34  
35   #undef  toupper
36   #define  toupper(c)     ((c) & ~0x20)   /* ASCII trick to convert case */
# Line 69 | Line 58 | struct {
58          {"", 0}
59   };
60                                          /* required values */
61 + int  year = 0;                                  /* year (optional) */
62   int  month, day;                                /* date */
63   double  hour;                                   /* time */
64   int  tsolar;                                    /* 0=standard, 1=solar */
# Line 78 | Line 68 | int  skytype = S_CLEAR;                                /* sky type */
68   int  dosun = 1;
69   double  zenithbr = 0.0;
70   int     u_zenith = 0;                           /* -1=irradiance, 1=radiance */
71 < double  turbidity = 2.75;
71 > double  turbidity = 2.45;
72   double  gprefl = 0.2;
73                                          /* computed values */
74   double  sundir[3];
# Line 95 | Line 85 | void printsky(void);
85   void printdefaults(void);
86   void userror(char  *msg);
87   double normsc(void);
88 < void cvthour(char  *hs);
99 < void printhead(register int  ac, register char  **av);
88 > int cvthour(char  *hs);
89  
90  
91   int
92 < main(argc, argv)
93 < int  argc;
94 < char  *argv[];
92 > main(
93 >        int  argc,
94 >        char  *argv[]
95 > )
96   {
97 +        int  got_meridian = 0;
98          int  i;
99  
100          progname = argv[0];
# Line 124 | Line 115 | char  *argv[];
115                  day = atoi(argv[2]);
116                  if (day < 1 || day > 31)
117                          userror("bad day");
118 <                cvthour(argv[3]);
118 >                got_meridian = cvthour(argv[3]);
119          }
120          for (i = 4; i < argc; i++)
121                  if (argv[i][0] == '-' || argv[i][0] == '+')
# Line 133 | Line 124 | char  *argv[];
124                                  skytype = S_CLEAR;
125                                  dosun = argv[i][0] == '+';
126                                  break;
127 +                        case 'y':
128 +                                year = atoi(argv[++i]);
129 +                                break;
130                          case 'r':
131                          case 'R':
132                                  u_solar = argv[i][1]=='R' ? -1 : 1;
# Line 166 | Line 160 | char  *argv[];
160                                  s_longitude = atof(argv[++i]) * (PI/180);
161                                  break;
162                          case 'm':
163 +                                if (got_meridian) {
164 +                                        ++i;
165 +                                        break;          /* time overrides */
166 +                                }
167                                  s_meridian = atof(argv[++i]) * (PI/180);
168                                  break;
169                          default:
# Line 175 | Line 173 | char  *argv[];
173                  else
174                          userror("bad option");
175  
176 <        if (fabs(s_meridian-s_longitude) > 45*PI/180)
176 >        if (year && (year < 1950) | (year > 2050))
177                  fprintf(stderr,
178 <        "%s: warning: %.1f hours btwn. standard meridian and longitude\n",
178 >                        "%s: warning - year should be in range 1950-2050\n",
179 >                        progname);
180 >        if (month && !tsolar && fabs(s_meridian-s_longitude) > 45*PI/180)
181 >                fprintf(stderr,
182 >        "%s: warning - %.1f hours btwn. standard meridian and longitude\n",
183                          progname, (s_longitude-s_meridian)*12/PI);
184  
185 <        printhead(argc, argv);
185 >        fputs("# ", stdout);
186 >        printargs(argc, argv, stdout);
187  
188          computesky();
189          printsky();
# Line 195 | Line 198 | computesky(void)                       /* compute sky parameters */
198          double  normfactor;
199                                          /* compute solar direction */
200          if (month) {                    /* from date and time */
201 <                int  jd;
199 <                double  sd, st;
201 >                double  sd, st = hour;
202  
203 <                jd = jdate(month, day);         /* Julian date */
204 <                sd = sdec(jd);                  /* solar declination */
205 <                if (tsolar)                     /* solar time */
206 <                        st = hour;
207 <                else
208 <                        st = hour + stadj(jd);
203 >                if (year) {                     /* Michalsky algorithm? */
204 >                        double  mjd = mjdate(year, month, day, hour);
205 >                        if (tsolar)
206 >                                sd = msdec(mjd, NULL);
207 >                        else
208 >                                sd = msdec(mjd, &st);
209 >                } else {
210 >                        int  jd = jdate(month, day);    /* Julian date */
211 >                        sd = sdec(jd);                  /* solar declination */
212 >                        if (!tsolar)                    /* get solar time? */
213 >                                st = hour + stadj(jd);
214 >                }
215                  altitude = salt(sd, st);
216                  azimuth = sazi(sd, st);
217                  printf("# Local solar time: %.2f\n", st);
# Line 281 | Line 289 | printsky(void)                 /* print out sky */
289          if (dosun) {
290                  printf("\nvoid light solar\n");
291                  printf("0\n0\n");
292 <                printf("3 %.2e %.2e %.2e\n", solarbr, solarbr, solarbr);
292 >                printf("3 %.3e %.3e %.3e\n", solarbr, solarbr, solarbr);
293                  printf("\nsolar source sun\n");
294                  printf("0\n0\n");
295                  printf("4 %f %f %f 0.5\n", sundir[0], sundir[1], sundir[2]);
# Line 291 | Line 299 | printsky(void)                 /* print out sky */
299          printf("2 skybr skybright.cal\n");
300          printf("0\n");
301          if (overcast)
302 <                printf("3 %d %.2e %.2e\n", skytype, zenithbr, groundbr);
302 >                printf("3 %d %.3e %.3e\n", skytype, zenithbr, groundbr);
303          else
304 <                printf("7 %d %.2e %.2e %.2e %f %f %f\n",
304 >                printf("7 %d %.3e %.3e %.3e %f %f %f\n",
305                                  skytype, zenithbr, groundbr, F2,
306                                  sundir[0], sundir[1], sundir[2]);
307   }
# Line 356 | Line 364 | normsc(void)                   /* compute normalization factor (E0*F2/
364                                  /* intermediate sky approx. */
365                  {3.5556, -2.7152, -1.3081, 1.0660, 0.60227},
366          };
367 <        register double  *nf;
367 >        double  *nf;
368          double  x, nsc;
369 <        register int  i;
369 >        int  i;
370                                          /* polynomial approximation */
371          nf = nfc[skytype==S_INTER];
372          x = (altitude - PI/4.0)/(PI/4.0);
# Line 370 | Line 378 | normsc(void)                   /* compute normalization factor (E0*F2/
378   }
379  
380  
381 < void
381 > int
382   cvthour(                        /* convert hour string */
383          char  *hs
384   )
385   {
386 <        register char  *cp = hs;
387 <        register int    i, j;
386 >        char  *cp = hs;
387 >        int     i, j;
388  
389          if ( (tsolar = *cp == '+') ) cp++;              /* solar time? */
390          while (isdigit(*cp)) cp++;
# Line 388 | Line 396 | cvthour(                       /* convert hour string */
396          }
397          while (isdigit(*cp)) cp++;
398          if (!*cp)
399 <                return;
399 >                return(0);
400          if (tsolar || !isalpha(*cp)) {
401                  fprintf(stderr, "%s: bad time format: %s\n", progname, hs);
402                  exit(1);
# Line 400 | Line 408 | cvthour(                       /* convert hour string */
408                                  break;
409                  if (!cp[j] && !tzone[i].zname[j]) {
410                          s_meridian = tzone[i].zmer * (PI/180);
411 <                        return;
411 >                        return(1);
412                  }
413          } while (tzone[i++].zname[0]);
414  
# Line 410 | Line 418 | cvthour(                       /* convert hour string */
418                  fprintf(stderr, " %s", tzone[i].zname);
419          putc('\n', stderr);
420          exit(1);
413 }
414
415
416 void
417 printhead(              /* print command header */
418        register int  ac,
419        register char  **av
420 )
421 {
422        putchar('#');
423        while (ac--) {
424                putchar(' ');
425                fputs(*av++, stdout);
426        }
427        putchar('\n');
421   }

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines