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

Comparing ray/src/gen/genssky.c (file contents):
Revision 2.4 by greg, Thu Aug 8 02:00:48 2024 UTC vs.
Revision 2.7 by greg, Thu Apr 10 23:30:58 2025 UTC

# Line 1 | Line 1
1 + #include "color.h"
2   #ifndef lint
3 < static const char RCSid[] = "$Id$";
3 > static const char RCSid[] =
4 >    "$Id$";
5   #endif
6   /* Main function for generating spectral sky */
7   /* Cloudy sky computed as weight average of clear and cie overcast sky */
# Line 160 | Line 162 | static inline double wmean(const double a, const doubl
162    return (a * x + b * y) / (a + b);
163   }
164  
165 < static double get_zenith_brightness(const double sundir[3]) {
165 > static double get_overcast_zenith_brightness(const double sundir[3]) {
166    double zenithbr;
167    if (sundir[2] < 0) {
168      zenithbr = 0;
# Line 190 | Line 192 | static void write_header(const int argc, char **argv,
192        cloud_cover, grefl, res);
193   }
194  
195 < static void write_rad(const double *sun_radiance, const FVECT sundir,
196 <                      const char *ddir, const char *skyfile) {
195 > static void write_rad(const double *sun_radiance, const double intensity,
196 >                      const FVECT sundir, const char *ddir,
197 >                      const char *skyfile) {
198    if (sundir[2] > 0) {
199      printf("void spectrum sunrad\n0\n0\n22 380 780 ");
197    /* Normalize to one */
198    double sum = 0.0;
200      int i;
201      for (i = 0; i < NSSAMP; ++i) {
202 <      sum += sun_radiance[i];
202 >      printf("%.3f ", sun_radiance[i]);
203      }
203    double mean = sum / NSSAMP;
204    for (i = 0; i < NSSAMP; ++i) {
205      printf("%.3f ", sun_radiance[i] / mean);
206    }
207    double intensity = mean * WVLSPAN;
204      printf("\n\nsunrad light solar\n0\n0\n3 %.1f %.1f %.1f\n\n", intensity,
205             intensity, intensity);
206      printf("solar source sun\n0\n0\n4 %f %f %f 0.533\n\n", sundir[0], sundir[1],
207             sundir[2]);
208    }
209 <  printf("void specpict skyfunc\n5 noop %s . 'atan2(Dy,Dx)/PI+1' "
210 <         "'acos(Dz)/PI'\n0\n0\n\n",
209 >  printf("void specpict skyfunc\n5 noop %s . 'Atan2(Dy,Dx)/PI+1' "
210 >         "'1-Acos(Dz)/PI'\n0\n0\n\n",
211           skyfile);
212   }
213  
# Line 231 | Line 227 | int gen_spect_sky(DATARRAY *tau_clear, DATARRAY *scat_
227                    DATARRAY *scat1m_clear, DATARRAY *irrad_clear,
228                    const double cloud_cover, const FVECT sundir,
229                    const double grefl, const int res, const char *outname,
230 <                  const char *ddir) {
230 >                  const char *ddir, const double dirnorm, const double difhor) {
231    char skyfile[PATH_MAX];
236  char grndfile[PATH_MAX];
232    if (!snprintf(skyfile, sizeof(skyfile), "%s%c%s_sky.hsr", ddir, DIRSEP,
233                  outname)) {
234      fprintf(stderr, "Error setting sky file name\n");
# Line 250 | Line 245 | int gen_spect_sky(DATARRAY *tau_clear, DATARRAY *scat_
245    FVECT view_point = {0, 0, ER + 10};
246    const double radius = VLEN(view_point);
247    const double sun_ct = fdot(view_point, sundir) / radius;
248 +
249 +  double overcast_zenithbr = get_overcast_zenith_brightness(sundir);
250 +  double overcast_grndbr = overcast_zenithbr * GNORM;
251 +
252 +  double dif_ratio = 1;
253 +  if (difhor > 0) {
254 +    DATARRAY *indirect_irradiance_clear = get_indirect_irradiance(irrad_clear, radius, sun_ct);
255 +    double overcast_ghi = overcast_zenithbr * 7.0 * PI / 9.0;
256 +    double diffuse_irradiance = 0;
257 +    int l;
258 +    for (l = 0; l < NSSAMP; ++l) {
259 +      diffuse_irradiance += indirect_irradiance_clear->arr.d[l] * 20;  /* 20nm interval */
260 +    }
261 +    free(indirect_irradiance_clear);
262 +    diffuse_irradiance = wmean2(diffuse_irradiance, overcast_ghi, cloud_cover);
263 +    if (diffuse_irradiance > 0) {
264 +        dif_ratio = difhor / WHTEFFICACY / diffuse_irradiance / 1.15;       /* fudge */
265 +    }
266 +  }
267    int i, j, k;
268    for (j = 0; j < yres; ++j) {
269      for (i = 0; i < xres; ++i) {
# Line 282 | Line 296 | int gen_spect_sky(DATARRAY *tau_clear, DATARRAY *scat_
296        }
297  
298        if (cloud_cover > 0) {
299 <        double zenithbr = get_zenith_brightness(sundir);
286 <        double grndbr = zenithbr * GNORM;
287 <        double skybr = get_overcast_brightness(rdir[2], zenithbr);
299 >        double skybr = get_overcast_brightness(rdir[2], overcast_zenithbr);
300          if (rdir[2] < 0) {
301            for (k = 0; k < NSSAMP; ++k) {
302 <            radiance[k] = wmean2(radiance[k], grndbr * D6415[k], cloud_cover);
302 >            radiance[k] = wmean2(radiance[k], overcast_grndbr * D6415[k], cloud_cover);
303            }
304          } else {
305            for (k = 0; k < NSSAMP; ++k) {
# Line 296 | Line 308 | int gen_spect_sky(DATARRAY *tau_clear, DATARRAY *scat_
308          }
309        }
310  
311 <      scolor2scolr(sky_sclr, radiance, 20);
311 >      for (k = 0; k < NSSAMP; ++k) {
312 >        radiance[k] *= dif_ratio;
313 >      }
314 >
315 >      scolor2scolr(sky_sclr, radiance, NSSAMP);
316        putbinary(sky_sclr, LSCOLR, 1, skyfp);
317      }
318    }
# Line 307 | Line 323 | int gen_spect_sky(DATARRAY *tau_clear, DATARRAY *scat_
323    get_solar_radiance(tau_clear, scat_clear, scat1m_clear, sundir, radius,
324                       sun_ct, sun_radiance);
325    if (cloud_cover > 0) {
326 <    double zenithbr = get_zenith_brightness(sundir);
311 <    double skybr = get_overcast_brightness(sundir[2], zenithbr);
326 >    double skybr = get_overcast_brightness(sundir[2], overcast_zenithbr);
327      int i;
328      for (i = 0; i < NSSAMP; ++i) {
329        sun_radiance[i] =
# Line 316 | Line 331 | int gen_spect_sky(DATARRAY *tau_clear, DATARRAY *scat_
331      }
332    }
333  
334 <  write_rad(sun_radiance, sundir, ddir, skyfile);
334 >  /* Normalize */
335 >  double sum = 0.0;
336 >  for (i = 0; i < NSSAMP; ++i) {
337 >    sum += sun_radiance[i];
338 >  }
339 >  double mean = sum / NSSAMP;
340 >  for (i = 0; i < NSSAMP; ++i) {
341 >    sun_radiance[i] /= mean;
342 >  }
343 >  double intensity = mean * WVLSPAN;
344 >  if (dirnorm > 0) {
345 >    intensity = dirnorm / SOLOMG / WHTEFFICACY;
346 >  }
347 >
348 >  write_rad(sun_radiance, intensity, sundir, ddir, skyfile);
349    return 1;
350   }
351  
# Line 419 | Line 448 | int main(int argc, char *argv[]) {
448    char lstag[3];
449    char *ddir = ".";
450    int i;
451 +  double dirnorm = 0; /* direct normal illuminance */
452 +  double difhor = 0;  /* diffuse horizontal illuminance */
453  
454    if (argc == 2 && !strcmp(argv[1], "-defaults")) {
455      printf("-i %d\t\t\t\t#scattering order\n", sorder);
# Line 434 | Line 465 | int main(int argc, char *argv[]) {
465    if (argc < 4) {
466      fprintf(stderr,
467              "Usage: %s month day hour -y year -a lat -o lon -m tz -d aod -r "
468 <            "res -n nproc -c ccover -l mie -g grefl -f outpath\n",
468 >            "res -n nproc -c ccover -l mie -L dirnorm_illum difhor_illum "
469 >            "-g grefl -f outpath\n",
470              argv[0]);
471      return 0;
472    }
# Line 494 | Line 526 | int main(int argc, char *argv[]) {
526        case 'o':
527          s_longitude = atof(argv[++i]) * (PI / 180.0);
528          break;
529 +      case 'L':
530 +        dirnorm = atof(argv[++i]);
531 +        difhor = atof(argv[++i]);
532 +        break;
533        case 'p':
534          ddir = argv[++i];
535          break;
# Line 563 | Line 599 | int main(int argc, char *argv[]) {
599    write_header(argc, argv, ccover, grefl, res);
600  
601    if (!gen_spect_sky(tau_clear_dp, scat_clear_dp, scat1m_clear_dp,
602 <                     irrad_clear_dp, ccover, sundir, grefl, res, outname,
603 <                     ddir)) {
602 >                     irrad_clear_dp, ccover, sundir, grefl, res, outname, ddir,
603 >                     dirnorm, difhor)) {
604      fprintf(stderr, "gen_spect_sky failed\n");
605      exit(1);
606    }

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines