| 89 |
|
#include "color.h" |
| 90 |
|
#include "sun.h" |
| 91 |
|
|
| 92 |
< |
char *progname; /* Program name */ |
| 93 |
< |
char errmsg[128]; /* Error message buffer */ |
| 92 |
> |
char *progname; /* Program name */ |
| 93 |
|
const double DC_SolarConstantE = 1367.0; /* Solar constant W/m^2 */ |
| 94 |
|
const double DC_SolarConstantL = 127.5; /* Solar constant klux */ |
| 95 |
|
|
| 134 |
|
/* Radiuans into degrees */ |
| 135 |
|
#define RadToDeg(rad) ((rad)*(180./PI)) |
| 136 |
|
|
| 138 |
– |
|
| 137 |
|
/* Perez sky model coefficients */ |
| 138 |
|
|
| 139 |
|
/* Reference: Perez, R., R. Seals, and J. Michalsky, 1993. "All- */ |
| 310 |
|
double rotation = 0; /* site rotation (degrees) */ |
| 311 |
|
double elevation; /* site elevation (meters) */ |
| 312 |
|
int leap_day = 0; /* add leap day? */ |
| 313 |
+ |
int sun_hours_only = 0; /* only output sun hours? */ |
| 314 |
|
int dir_is_horiz; /* direct is meas. on horizontal? */ |
| 315 |
|
FILE *sunsfp = NULL; /* output file for individual suns */ |
| 316 |
|
FILE *modsfp = NULL; /* modifier output file */ |
| 405 |
|
case 's': /* sky only (no direct) */ |
| 406 |
|
suncolor[0] = suncolor[1] = suncolor[2] = 0; |
| 407 |
|
break; |
| 408 |
+ |
case 'u': /* solar hours only */ |
| 409 |
+ |
sun_hours_only = 1; |
| 410 |
+ |
break; |
| 411 |
|
case 'r': /* rotate distribution */ |
| 412 |
|
if (argv[i][2] && argv[i][2] != 'z') |
| 413 |
|
goto userr; |
| 500 |
|
/* process each time step in tape */ |
| 501 |
|
while (scanf("%d %d %lf %lf %lf\n", &mo, &da, &hr, &dir, &dif) == 5) { |
| 502 |
|
double sda, sta; |
| 503 |
< |
|
| 502 |
< |
mtx_offset = 3*nskypatch*nstored; |
| 503 |
< |
nstored += !avgSky | !nstored; |
| 504 |
< |
/* make space for next row */ |
| 505 |
< |
if (nstored > tstorage) { |
| 506 |
< |
tstorage += (tstorage>>1) + nstored + 7; |
| 507 |
< |
mtx_data = resize_dmatrix(mtx_data, tstorage, nskypatch); |
| 508 |
< |
} |
| 509 |
< |
ntsteps++; /* keep count of time steps */ |
| 503 |
> |
int sun_in_sky; |
| 504 |
|
/* compute solar position */ |
| 505 |
|
if ((mo == 2) & (da == 29)) { |
| 506 |
|
julian_date = 60; |
| 510 |
|
sda = sdec(julian_date); |
| 511 |
|
sta = stadj(julian_date); |
| 512 |
|
altitude = salt(sda, hr+sta); |
| 513 |
+ |
sun_in_sky = (altitude > -DegToRad(SUN_ANG_DEG/2.)); |
| 514 |
+ |
if (sun_hours_only && !sun_in_sky) |
| 515 |
+ |
continue; /* skipping nighttime points */ |
| 516 |
|
azimuth = sazi(sda, hr+sta) + PI - DegToRad(rotation); |
| 517 |
|
|
| 518 |
+ |
mtx_offset = 3*nskypatch*nstored; |
| 519 |
+ |
nstored += !avgSky | !nstored; |
| 520 |
+ |
/* make space for next row */ |
| 521 |
+ |
if (nstored > tstorage) { |
| 522 |
+ |
tstorage += (tstorage>>1) + nstored + 7; |
| 523 |
+ |
mtx_data = resize_dmatrix(mtx_data, tstorage, nskypatch); |
| 524 |
+ |
} |
| 525 |
+ |
ntsteps++; /* keep count of time steps */ |
| 526 |
+ |
|
| 527 |
|
if (dir+dif <= 1e-4) { /* effectively nighttime? */ |
| 528 |
|
if (!avgSky | !mtx_offset) |
| 529 |
|
memset(mtx_data+mtx_offset, 0, |
| 530 |
|
sizeof(float)*3*nskypatch); |
| 531 |
< |
if (sunsfp) /* output black sun */ |
| 531 |
> |
/* output black sun? */ |
| 532 |
> |
if (sunsfp && sun_in_sky) |
| 533 |
|
OutputSun(solar_minute(julian_date,hr), 0, |
| 534 |
|
sunsfp, modsfp); |
| 535 |
|
continue; |
| 536 |
|
} |
| 537 |
+ |
if (!sun_in_sky && dir > (input==1 ? 20. : 20.*WHTEFFICACY)) |
| 538 |
+ |
fprintf(stderr, |
| 539 |
+ |
"%s: warning - unusually bright at %.1f on %d-%d\n", |
| 540 |
+ |
progname, hr, mo, da); |
| 541 |
|
/* convert measured values */ |
| 542 |
< |
if (dir_is_horiz && altitude > 0.) |
| 542 |
> |
if (dir_is_horiz && altitude > FTINY) |
| 543 |
|
dir /= sin(altitude); |
| 544 |
|
if (input == 1) { |
| 545 |
|
dir_irrad = dir; |
| 550 |
|
} |
| 551 |
|
/* compute sky patch values */ |
| 552 |
|
ComputeSky(mtx_data+mtx_offset); |
| 553 |
< |
|
| 554 |
< |
if (sunsfp) /* output sun if requested */ |
| 553 |
> |
/* output sun if requested */ |
| 554 |
> |
if (sunsfp && sun_in_sky) |
| 555 |
|
OutputSun(solar_minute(julian_date,hr), 1, |
| 556 |
|
sunsfp, modsfp); |
| 557 |
|
|
| 653 |
|
fprintf(stderr, "%s: done.\n", progname); |
| 654 |
|
exit(0); |
| 655 |
|
userr: |
| 656 |
< |
fprintf(stderr, "Usage: %s [-v][-h][-A][-d|-s|-n][-D file][-r deg][-m N][-g r g b][-c r g b][-o{f|d}][-O{0|1}] [tape.wea]\n", |
| 656 |
> |
fprintf(stderr, "Usage: %s [-v][-h][-A][-d|-s|-n][-u][-D file [-M modfile]][-r deg][-m N][-g r g b][-c r g b][-o{f|d}][-O{0|1}] [tape.wea]\n", |
| 657 |
|
progname); |
| 658 |
|
exit(1); |
| 659 |
|
fmterr: |
| 825 |
|
FVECT sv; |
| 826 |
|
|
| 827 |
|
srad = DegToRad(SUN_ANG_DEG/2.); |
| 817 |
– |
|
| 818 |
– |
if (altitude < -srad) /* well below horizon? */ |
| 819 |
– |
return; |
| 820 |
– |
|
| 828 |
|
srad = goodsun ? dir_illum/(WHTEFFICACY * PI*srad*srad) : 0; |
| 829 |
|
vector(sv, altitude, azimuth); |
| 830 |
|
fprintf(fp, "\nvoid light solar%d\n0\n0\n", id); |
