| 308 |
|
if (ccover > 0) { |
| 309 |
|
double zenithbr = get_zenith_brightness(sundir); |
| 310 |
|
double skybr = get_overcast_brightness(sundir[2], zenithbr); |
| 311 |
< |
for (int l = 0; l < NSSAMP; ++l) { |
| 311 |
> |
int l; |
| 312 |
> |
for (l = 0; l < NSSAMP; ++l) { |
| 313 |
|
sun_radiance[l] = |
| 314 |
|
wmean2(sun_radiance[l], D6415[l] * skybr / WVLSPAN, ccover); |
| 315 |
|
} |
| 321 |
|
/* add to nearest patch radiances */ |
| 322 |
|
for (i = nsuns; i--;) { |
| 323 |
|
float *pdest = parr + NSSAMP * near_patch[i]; |
| 324 |
< |
for (int k = 0; k < NSSAMP; k++) { |
| 324 |
> |
int k; |
| 325 |
> |
for (k = 0; k < NSSAMP; k++) { |
| 326 |
|
*pdest++ = sun_radiance[k] * wta[i] / wtot; |
| 327 |
|
} |
| 328 |
|
} |
| 333 |
|
int i; |
| 334 |
|
double mu_sky; /* Sun-sky point azimuthal angle */ |
| 335 |
|
double sspa; /* Sun-sky point angle */ |
| 334 |
– |
double zsa; /* Zenithal sun angle */ |
| 336 |
|
FVECT view_point = {0, 0, ER}; |
| 337 |
|
for (i = 1; i < nskypatch; i++) { |
| 338 |
|
FVECT rdir_sky; |
| 339 |
+ |
int k; |
| 340 |
|
vectorize(rh_palt[i], rh_pazi[i], rdir_sky); |
| 341 |
|
mu_sky = fdot(view_point, rdir_sky) / ER; |
| 342 |
|
sspa = fdot(rdir_sky, sundir); |
| 343 |
|
SCOLOR sky_radiance = {0}; |
| 344 |
|
|
| 345 |
|
get_sky_radiance(scat, scat1m, ER, mu_sky, sun_ct, sspa, sky_radiance); |
| 346 |
< |
for (int k = 0; k < NSSAMP; ++k) { |
| 346 |
> |
for (k = 0; k < NSSAMP; ++k) { |
| 347 |
|
sky_radiance[k] *= WVLSPAN; |
| 348 |
|
} |
| 349 |
|
|
| 351 |
|
double zenithbr = get_zenith_brightness(sundir); |
| 352 |
|
double grndbr = zenithbr * GNORM; |
| 353 |
|
double skybr = get_overcast_brightness(rdir_sky[2], zenithbr); |
| 354 |
< |
for (int k = 0; k < NSSAMP; ++k) { |
| 354 |
> |
int k; |
| 355 |
> |
for (k = 0; k < NSSAMP; ++k) { |
| 356 |
|
sky_radiance[k] = wmean2(sky_radiance[k], skybr * D6415[k], ccover); |
| 357 |
|
} |
| 358 |
|
} |
| 359 |
|
|
| 360 |
< |
for (int k = 0; k < NSSAMP; ++k) { |
| 360 |
> |
for (k = 0; k < NSSAMP; ++k) { |
| 361 |
|
parr[NSSAMP * i + k] = sky_radiance[k]; |
| 362 |
|
} |
| 363 |
|
} |
| 382 |
|
const FVECT rdir_grnd = {0, 0, -1}; |
| 383 |
|
const double mu_grnd = fdot(view_point, rdir_grnd) / radius; |
| 384 |
|
const double nu_grnd = fdot(rdir_grnd, sundir); |
| 385 |
+ |
int j; |
| 386 |
|
|
| 387 |
|
/* Calculate sun zenith angle (don't let it dip below horizon) */ |
| 388 |
|
/* Also limit minimum angle to keep circumsolar off zenith */ |
| 396 |
|
/* Compute ground radiance (include solar contribution if any) */ |
| 397 |
|
get_ground_radiance(tau, scat, scat1m, irrad, view_point, rdir_grnd, radius, |
| 398 |
|
mu_grnd, sun_ct, nu_grnd, grefl, sundir, parr); |
| 399 |
< |
for (int j = 0; j < NSSAMP; j++) { |
| 399 |
> |
for (j = 0; j < NSSAMP; j++) { |
| 400 |
|
parr[j] *= WVLSPAN; |
| 401 |
|
} |
| 402 |
|
calc_sky_patch_radiance(scat, scat1m, ccover, parr); |
| 542 |
|
|
| 543 |
|
while (scanf("%d %d %lf %lf %lf %lf %lf\n", &mo, &da, &hr, &dni, &dhi, &aod, |
| 544 |
|
&cc) == 7) { |
| 545 |
+ |
if (aod == 0.0) { |
| 546 |
+ |
aod = AOD0_CA; |
| 547 |
+ |
fprintf(stderr, "aod is zero, using default value %.3f\n", AOD0_CA); |
| 548 |
+ |
} |
| 549 |
|
double sda, sta; |
| 550 |
|
int sun_in_sky; |
| 551 |
|
/* compute solar position */ |
| 569 |
|
|
| 570 |
|
mtx_offset = NSSAMP * nskypatch * nstored; |
| 571 |
|
nstored += 1; |
| 572 |
+ |
printf("mtx_offset = %d nstored = %d nskypatch = %d\n", mtx_offset, nstored, |
| 573 |
+ |
nskypatch); |
| 574 |
|
/* make space for next row */ |
| 575 |
|
if (nstored > tstorage) { |
| 576 |
+ |
printf("make space for next row nstored = %d tstorage = %d\n", nstored, |
| 577 |
+ |
tstorage); |
| 578 |
|
tstorage += (tstorage >> 1) + nstored + 7; |
| 579 |
|
mtx_data = resize_dmatrix(mtx_data, tstorage, nskypatch); |
| 580 |
|
} |
| 622 |
|
if (!sky_only) |
| 623 |
|
add_direct(tau_clear_dp, scat_clear_dp, scat1m_clear_dp, irrad_clear_dp, |
| 624 |
|
cc, mtx_data + mtx_offset); |
| 613 |
– |
/* update cumulative sky? */ |
| 614 |
– |
for (i = NSSAMP * nskypatch * (ntsteps > 1); i--;) |
| 615 |
– |
mtx_data[i] += mtx_data[mtx_offset + i]; |
| 625 |
|
/* monthly reporting */ |
| 626 |
|
if (verbose && mo != last_monthly) |
| 627 |
|
fprintf(stderr, "%s: stepping through month %d...\n", progname, |
| 671 |
|
switch (outfmt) { |
| 672 |
|
case 'a': |
| 673 |
|
for (j = 0; j < nstored; j++) { |
| 674 |
< |
for (int k = 0; k < NSSAMP; k++) { |
| 675 |
< |
printf("%.3g \n", mtx_data[mtx_offset + k]); |
| 674 |
> |
int k; |
| 675 |
> |
for (k = 0; k < NSSAMP; k++) { |
| 676 |
> |
printf("%.3g ", mtx_data[mtx_offset + k]); |
| 677 |
|
} |
| 678 |
|
printf("\n"); |
| 679 |
|
mtx_offset += NSSAMP * nskypatch; |