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

Comparing src/gen/gendaylit.c (file contents):
Revision 2.12 by greg, Fri Aug 9 16:44:19 2013 UTC vs.
Revision 2.15 by greg, Mon Nov 18 18:07:16 2013 UTC

# Line 6 | Line 6
6   *                              1 Avenue Eugene Freyssinet, Saint-Quentin-Yvelines, France
7   */
8  
9 + #define  _USE_MATH_DEFINES
10   #include  <stdio.h>
11   #include  <string.h>
12   #include  <math.h>
# Line 55 | Line 56 | float defangle_phi[] = {
56  
57  
58  
59 < /* Perez sky parametrization : epsilon and delta calculations from the direct and diffuse irradiances */
59 > /* Perez sky parametrization: epsilon and delta calculations from the direct and diffuse irradiances */
60   double sky_brightness();
61   double sky_clearness();
62  
# Line 88 | Line 89 | void printdefaults();
89   void check_sun_position();
90   void computesky();
91   void printhead(int ac, char** av);
92 < void userror(char* msg);
92 > void usage_error(char* msg);
93   void printsky();
94  
95   FILE * frlibopen(char* fname);
# Line 97 | Line 98 | FILE * frlibopen(char* fname);
98   double  get_eccentricity();
99   double  air_mass();
100  
101 < double  solar_sunset();
102 < double  solar_sunrise();
101 > double  solar_sunset(int month, int day);
102 > double  solar_sunrise(int month, int day);
103   double  stadj();
104   int     jdate(int month, int day);
105  
106  
106 –
107   /* sun calculation constants */
108   extern double   s_latitude;
109   extern double   s_longitude;
# Line 111 | Line 111 | extern double   s_meridian;
111  
112   const double    AU = 149597890E3;
113   const double    solar_constant_e = 1367;    /* solar constant W/m^2 */
114 < const double    solar_constant_l = 127.5;   /* solar constant klux */
114 > const double    solar_constant_l = 127500;   /* solar constant lux */
115  
116   const double    half_sun_angle = 0.2665;
117   const double    half_direct_angle = 2.85;
118  
119 < const double    skyclearinf = 1.0;      /* limitations for the variation of the Perez parameters */
120 < const double    skyclearsup = 12.1;
119 > const double    skyclearinf = 1.0;          /* limitations for the variation of the Perez parameters */
120 > const double    skyclearsup = 12.01;
121   const double    skybriginf = 0.01;
122   const double    skybrigsup = 0.6;
123  
# Line 168 | Line 168 | float  timeinterval = 0;
168   char    *progname;
169   char    errmsg[128];
170  
171 + double  st;
172  
173  
173 –
174   int main(int argc, char** argv)
175   {
176          int  i;
# Line 181 | Line 181 | int main(int argc, char** argv)
181                  return 0;
182          }
183          if (argc < 4)
184 <                userror("arg count");
184 >                usage_error("arg count");
185          if (!strcmp(argv[1], "-ang")) {
186                  altitude = atof(argv[2]) * (M_PI/180);
187                  azimuth = atof(argv[3]) * (M_PI/180);
# Line 189 | Line 189 | int main(int argc, char** argv)
189          } else {
190                  month = atoi(argv[1]);
191                  if (month < 1 || month > 12)
192 <                        userror("bad month");
192 >                        usage_error("bad month");
193                  day = atoi(argv[2]);
194                  if (day < 1 || day > 31)
195 <                        userror("bad day");
195 >                        usage_error("bad day");
196                  hour = atof(argv[3]);
197                  if (hour < 0 || hour >= 24)
198 <                        userror("bad hour");
198 >                        usage_error("bad hour");
199                  tsolar = argv[3][0] == '+';
200          }
201          for (i = 4; i < argc; i++)
# Line 236 | Line 236 | int main(int argc, char** argv)
236                                  break;
237                          
238                          case 'O':
239 <                                output = atof(argv[++i]);       /*define the unit of the output of the program :
240 <                                                                sky and sun luminance/radiance (0==W visible, 1==W solar radiation, 2==lm) */
239 >                                output = atof(argv[++i]);       /*define the unit of the output of the program:
240 >                                                                sky and sun luminance/radiance
241 >                                                                (0==W visible, 1==W solar radiation, 2==lm) */
242                                  break;
243                                  
244                          case 'P':
# Line 269 | Line 270 | int main(int argc, char** argv)
270                                  globalirradiance = atof(argv[++i]);
271                                  break;
272                          
272 –                        /*
273 –                        case 'l':
274 –                                sunaltitude_border = atof(argv[++i]);
275 –                                break;
276 –                        */
277 –                                
273                          case 'i':
274                                  timeinterval = atof(argv[++i]);
275                                  break;
# Line 282 | Line 277 | int main(int argc, char** argv)
277                          
278                          default:
279                                  sprintf(errmsg, "unknown option: %s", argv[i]);
280 <                                userror(errmsg);
280 >                                usage_error(errmsg);
281                          }
282                  else
283 <                        userror("bad option");
283 >                        usage_error("bad option");
284  
285 <        if (fabs(s_meridian-s_longitude) > 30*M_PI/180)
286 <                fprintf(stderr,
292 <                    "%s: warning: %.1f hours btwn. standard meridian and longitude\n",
285 >        if (month && !tsolar && fabs(s_meridian-s_longitude) > 45*M_PI/180)
286 >                fprintf(stderr,"%s: warning: %.1f hours btwn. standard meridian and longitude\n",
287                      progname, (s_longitude-s_meridian)*12/M_PI);
288  
289  
290          /* dynamic memory allocation for the pointers */
297 –        
291          if ( (c_perez = calloc(5, sizeof(double))) == NULL )
292 <        {
300 <                fprintf(stderr,"Out of memory error in function main");
301 <                return 1;
302 <        }
292 >        { fprintf(stderr,"Out of memory error in function main"); return 1; }
293  
294 +        
295          printhead(argc, argv);
296          computesky();
306 –        
307 –        if(*(c_perez+1)>0)
308 –        {
309 –          fprintf(stderr, "Warning: positive Perez parameter B (= %lf), printing error sky\n",*(c_perez+1));      
310 –          print_error_sky();
311 –          exit(1);
312 –        }
313 –        
297          printsky();
298          return 0;
299 +
300   }
301  
302  
303  
304  
305  
322 –
323 –
324 –
306   void computesky()
307   {
308  
# Line 338 | Line 319 | void computesky()
319                  
320          if (month) {                    /* from date and time */
321                  int  jd;
322 <                double  sd, st;
322 >                double  sd;
323  
324                  jd = jdate(month, day);         /* Julian date */
325                  sd = sdec(jd);                  /* solar declination */
# Line 346 | Line 327 | void computesky()
327                          st = hour;
328                  else
329                          st = hour + stadj(jd);
349 –                                
330                  
351 –                if(st<solar_sunrise(month,day) || st>solar_sunset(month,day)) {
352 –                print_error_sky();
353 –                exit(1);
354 –                }
355 –                
331                                          
332                  if(timeinterval) {
333                          
# Line 360 | Line 335 | void computesky()
335                          fprintf(stderr, "time interval negative\n");
336                          exit(1);
337                          }
338 <                                        
339 <                        if(fabs(solar_sunrise(month,day)-st)<timeinterval/60) {                
340 <                        
341 <                        fprintf(stderr, "Solar position corrected at %d %d %.3f\n",month,day,hour);
342 <                        st= (st+timeinterval/120+solar_sunrise(month,day))/2;
338 >                                                                        
339 >                        if(fabs(solar_sunrise(month,day)-st)<=timeinterval/120) {                      
340 >                         st= (st+timeinterval/120+solar_sunrise(month,day))/2;
341 >                         if(suppress_warnings==0)
342 >                         { fprintf(stderr, "Solar position corrected at time step %d %d %.3f\n",month,day,hour); }
343                          }
344                  
345 <                        if(fabs(solar_sunset(month,day)-st)<timeinterval/60) {
346 <                        fprintf(stderr, "Solar position corrected at %d %d %.3f\n",month,day,hour);
347 <                        st= (st-timeinterval/120+solar_sunset(month,day))/2;
345 >                        if(fabs(solar_sunset(month,day)-st)<timeinterval/120) {
346 >                         st= (st-timeinterval/120+solar_sunset(month,day))/2;
347 >                         if(suppress_warnings==0)
348 >                         { fprintf(stderr, "Solar position corrected at time step %d %d %.3f\n",month,day,hour); }
349                          }
350 +                        
351 +                        if((st<solar_sunrise(month,day)-timeinterval/120) || (st>solar_sunset(month,day)+timeinterval/120)) {
352 +                          if(suppress_warnings==0)
353 +                          { fprintf(stderr, "Warning: sun position too low, printing error sky at %d %d %.3f\n",month,day,hour); }
354 +                         altitude = salt(sd, st);
355 +                         azimuth = sazi(sd, st);
356 +                         print_error_sky();
357 +                         exit(0);
358 +                        }
359                  }
360 +                else
361                  
362 +                if(st<solar_sunrise(month,day) || st>solar_sunset(month,day)) {
363 +                        if(suppress_warnings==0)
364 +                        { fprintf(stderr, "Warning: sun altitude below zero at time step %i %i %.2f, printing error sky\n",month,day,hour); }
365 +                        altitude = salt(sd, st);
366 +                        azimuth = sazi(sd, st);
367 +                        print_error_sky();
368 +                        exit(0);
369 +                }
370                  
371                  altitude = salt(sd, st);
372                  azimuth = sazi(sd, st);
# Line 384 | Line 378 | void computesky()
378          
379  
380  
387 –        /* if loop for the -l option. W.Sprenger (01/2013) */
388 –        /*
389 –        if (altitude*180/M_PI < sunaltitude_border) {
390 –        
391 –        if (suppress_warnings==0) {
392 –            fprintf(stderr, "Warning: sun altitude (%.3f degrees) below the border (%.3f degrees)\n",altitude*180/M_PI,sunaltitude_border);
393 –                                   }
394 –        print_error_sky();
395 –        exit(1);
396 –        }
397 –        */
398 –
381                          
382          if (!cloudy && altitude > 87.*M_PI/180.) {
383                  
# Line 407 | Line 389 | void computesky()
389                  altitude = 87.*M_PI/180.;
390          }
391          
392 +        
393 +        
394          sundir[0] = -sin(azimuth)*cos(altitude);
395          sundir[1] = -cos(azimuth)*cos(altitude);
396          sundir[2] = sin(altitude);
# Line 414 | Line 398 | void computesky()
398                  
399          /* calculation for the new functions */
400          sunzenith = 90 - altitude*180/M_PI;
401 <        
418 <        
401 >                        
402  
403          /* compute the inputs for the calculation of the light distribution over the sky*/
404          if (input==0)           /* P */
# Line 439 | Line 422 | void computesky()
422                  check_irradiances();
423                  skybrightness = sky_brightness();
424                  skyclearness =  sky_clearness();
425 +                
426                  check_parametrization();
427 <
427 >                                                        
428                  if (output==0 || output==2)
429                          {
430                          diffuseilluminance = diffuseirradiance*glob_h_diffuse_effi_PEREZ();/*diffuse horizontal illuminance*/
# Line 464 | Line 448 | void computesky()
448                          if (altitude<=0)
449                          {
450                                  if (suppress_warnings==0)
451 <                                     fprintf(stderr, "Warning: solar zenith angle larger than 90 degrees; using zero irradiance to proceed\n");
451 >                                     fprintf(stderr, "Warning: sun altitude < 0, proceed with irradiance values of zero\n");
452                                  directirradiance = 0;
453                                  diffuseirradiance = 0;
454                          } else {
# Line 504 | Line 488 | void computesky()
488                          {
489                                  if (suppress_warnings==0)
490                                          fprintf(stderr, "Warning: global irradiance is higher than the time-dependent solar constant s0\n");
507 –                                
491                                  globalirradiance=erbs_s0*0.999;                
492                          }
493                          
# Line 517 | Line 500 | void computesky()
500                          directirradiance=globalirradiance-diffuseirradiance;
501                          
502                          printf("# erbs_s0, erbs_kt, irr_dir_h, irr_diff: %.3f %.3f %.3f %.3f\n", erbs_s0, erbs_kt, directirradiance, diffuseirradiance);
503 <                        printf("# WARNING: the -E option is only recommended for a rough estimation!");
503 >                        printf("# WARNING: the -E option is only recommended for a rough estimation!\n");
504                          
505                          directirradiance=directirradiance/sin(altitude);
506                                                                                                                          
# Line 540 | Line 523 | void computesky()
523                  
524                  
525          
526 <        else    {fprintf(stderr,"error at the input arguments"); exit(1);}
526 >        else    { fprintf(stderr,"error at the input arguments"); exit(1); }
527  
528  
529          
# Line 562 | Line 545 | void computesky()
545  
546          
547          
565 –        
566 –        
548          /*calculation of the modelled luminance */
549          for (j=0;j<145;j++)
550          {
# Line 609 | Line 590 | void computesky()
590          else
591                  solarradiance = directilluminance/(2*M_PI*(1-cos(half_sun_angle*M_PI/180)));
592          
612 –                        
593  
594  
595 < /* Compute the ground radiance */
596 < zenithbr=calc_rel_lum_perez(0.0,radians(sunzenith),radians(sunzenith),skyclearness,skybrightness,coeff_perez);
597 < zenithbr*=diffnormalization;
595 >        /* Compute the ground radiance */
596 >        zenithbr=calc_rel_lum_perez(0.0,radians(sunzenith),radians(sunzenith),skyclearness,skybrightness,coeff_perez);
597 >        zenithbr*=diffnormalization;
598          
599 < if (skyclearness==1)
599 >        if (skyclearness==1)
600          normfactor = 0.777778;
601                  
602 < if (skyclearness>=6)
602 >        if (skyclearness>=6)
603          {              
604          F2 = 0.274*(0.91 + 10.0*exp(-3.0*(M_PI/2.0-altitude)) + 0.45*sundir[2]*sundir[2]);
605          normfactor = normsc()/F2/M_PI;
606          }
607  
608 < if ( (skyclearness>1) && (skyclearness<6) )
608 >        if ( (skyclearness>1) && (skyclearness<6) )
609          {
610          S_INTER=1;
611          F2 = (2.739 + .9891*sin(.3119+2.6*altitude)) * exp(-(M_PI/2.0-altitude)*(.4441+1.48*altitude));
612          normfactor = normsc()/F2/M_PI;
613          }
614  
615 < groundbr = zenithbr*normfactor;
615 >        groundbr = zenithbr*normfactor;
616  
617 < if (dosun&&(skyclearness>1))
617 >        if (dosun&&(skyclearness>1))
618          groundbr += 6.8e-5/M_PI*solarradiance*sundir[2];                
619  
620 < groundbr *= gprefl;
620 >        groundbr *= gprefl;
621  
622  
623 +                
624 +        if(*(c_perez+1)>0)
625 +        {
626 +          if(suppress_warnings==0)
627 +                {  fprintf(stderr, "Warning: positive Perez parameter B (= %lf), printing error sky\n",*(c_perez+1));}  
628 +          print_error_sky();
629 +          exit(0);
630 +        }
631  
632 +
633   return;
634   }
635  
636  
637  
638  
650 – void print_error_sky()
651 – {
652 –        sundir[0] = -sin(azimuth)*cos(altitude);
653 –        sundir[1] = -cos(azimuth)*cos(altitude);
654 –        sundir[2] = sin(altitude);
655 –        
656 –        printf("\nvoid brightfunc skyfunc\n");
657 –        printf("2 skybright perezlum.cal\n");
658 –        printf("0\n");
659 –        printf("10 0.00 0.00  0.000 0.000 0.000 0.000 0.000  %f %f %f \n", sundir[0], sundir[1], sundir[2]);
660 – }
661 –        
639  
663 –
664 –
665 –
666 –
667 –
640   double solar_sunset(int month,int day)
641   {
642       float W;
# Line 674 | Line 646 | double solar_sunset(int month,int day)
646   }
647  
648  
649 +
650 +
651   double solar_sunrise(int month,int day)
652   {
653       float W;
# Line 685 | Line 659 | double solar_sunrise(int month,int day)
659  
660  
661  
662 + void printsky()
663 + {      
664 +        
665 +        printf("# Local solar time: %.2f\n", st);
666 +        printf("# Solar altitude and azimuth: %.1f %.1f\n", altitude*180/M_PI, azimuth*180/M_PI);
667  
668  
690 –
691 –
692 –
693 –
694 –
695 – void printsky()                 /* print out sky */
696 – {
669          if (dosun&&(skyclearness>1))
670          {              
671                  printf("\nvoid light solar\n");
# Line 722 | Line 694 | void printsky()                        /* print out sky */
694   }
695  
696  
697 +
698 + void print_error_sky()
699 + {
700 +
701 +
702 +        sundir[0] = -sin(azimuth)*cos(altitude);
703 +        sundir[1] = -cos(azimuth)*cos(altitude);
704 +        sundir[2] = sin(altitude);
705 +
706 +        printf("# Local solar time: %.2f\n", st);
707 +        printf("# Solar altitude and azimuth: %.1f %.1f\n", altitude*180/M_PI, azimuth*180/M_PI);
708 +
709 +        printf("\nvoid brightfunc skyfunc\n");
710 +        printf("2 skybright perezlum.cal\n");
711 +        printf("0\n");
712 +        printf("10 0.00 0.00  0.000 0.000 0.000 0.000 0.000  %f %f %f \n", sundir[0], sundir[1], sundir[2]);
713 + }
714 +        
715 +
716 +
717 +
718 +
719   void printdefaults()                    /* print default values */
720   {
721          printf("-g %f\t\t\t# Ground plane reflectance\n", gprefl);
# Line 735 | Line 729 | void printdefaults()                   /* print default values */
729   }
730  
731  
732 < void userror(char* msg)                 /* print usage error and quit */
732 >
733 >
734 > void usage_error(char* msg)                     /* print usage error and quit */
735   {
736          if (msg != NULL)
737                  fprintf(stderr, "%s: Use error - %s\n\n", progname, msg);
# Line 752 | Line 748 | void userror(char* msg)                        /* print usage error and qui
748          fprintf(stderr, "       -E global-horizontal-irradiance (W/m^2)\n\n");
749          fprintf(stderr, "       Output specification with option:\n");
750          fprintf(stderr, "       -O [0|1|2]  (0=output in W/m^2/sr visible, 1=output in W/m^2/sr solar, 2=output in candela/m^2), default is 0 \n");
751 <        fprintf(stderr, "       gendaylit version 2.3 (2013/08/08)  \n\n");
751 >        fprintf(stderr, "       gendaylit version 2.4 (2013/09/04)  \n\n");
752          exit(1);
753   }
754  
755  
756  
757 +
758   double normsc()           /* compute normalization factor (E0*F2/L0) */
759   {
760          static double  nfc[2][5] = {
# Line 781 | Line 778 | double normsc()                  /* compute normalization fac
778  
779  
780  
781 +
782 +
783   void printhead(int ac, char** av)               /* print command header */
784   {
785          putchar('#');
# Line 794 | Line 793 | void printhead(int ac, char** av)              /* print command he
793  
794  
795  
796 +
797 +
798   /* Perez models */
799  
800   /* Perez global horizontal luminous efficacy model */
# Line 865 | Line 866 | double glob_h_effi_PEREZ()
866  
867  
868  
868 –
869          for (i=1; i<=category_total_number; i++)
870          {
871                  if ( (skyclearness >= category_bounds[i]) && (skyclearness < category_bounds[i+1]) )
# Line 879 | Line 879 | double glob_h_effi_PEREZ()
879   }
880  
881  
882 +
883 +
884   /* global horizontal diffuse efficacy model, according to PEREZ */
885   double glob_h_diffuse_effi_PEREZ()
886   {
# Line 886 | Line 888 | double glob_h_diffuse_effi_PEREZ()
888          double  category_bounds[10], a[10], b[10], c[10], d[10];
889          int     category_total_number, category_number, i;
890  
889 –
890 –
891 –        
891          check_parametrization();
892          
893          
894   /*if ((skyclearness<skyclearinf || skyclearness>skyclearsup || skybrightness<skybriginf || skybrightness>skybrigsup) && suppress_warnings==0)
895 <     fprintf(stderr, "Warning: skyclearness or skybrightness out of range in function glob_h_diffuse_PEREZ \n"); */
895 > fprintf(stderr, "Warning: skyclearness or skybrightness out of range in function glob_h_diffuse_PEREZ \n"); */
896      
898 –
899 –
897   /* initialize category bounds (clearness index bounds) */
898  
899          category_total_number = 8;
# Line 952 | Line 949 | double glob_h_diffuse_effi_PEREZ()
949  
950  
951  
955 –
952          category_number = -1;
953          for (i=1; i<=category_total_number; i++)
954          {
# Line 962 | Line 958 | double glob_h_diffuse_effi_PEREZ()
958  
959          if (category_number == -1) {
960                  if (suppress_warnings==0)
961 <                    fprintf(stderr, "ERROR: Model parameters out of range, skyclearness = %lf \n", skyclearness);
961 >                    fprintf(stderr, "Warning: sky clearness (= %.3f) too high, printing error sky\n", skyclearness);
962                  print_error_sky();
963 <                exit(1);
963 >                exit(0);
964          }
965                  
966  
# Line 977 | Line 973 | double glob_h_diffuse_effi_PEREZ()
973  
974  
975  
976 +
977 +
978 +
979   /* direct normal efficacy model, according to PEREZ */
980  
981   double direct_n_effi_PEREZ()
# Line 987 | Line 986 | double         category_bounds[10], a[10], b[10], c[10], d[10
986   int     category_total_number, category_number, i;
987  
988  
989 < if ((skyclearness<skyclearinf || skyclearness>skyclearsup || skybrightness<skybriginf || skybrightness>skybrigsup) && suppress_warnings==0)
990 <   fprintf(stderr, "Warning: skyclearness or skybrightness out of range in function direct_n_effi_PEREZ \n");
989 > /*if ((skyclearness<skyclearinf || skyclearness>skyclearsup || skybrightness<skybriginf || skybrightness>skybrigsup) && suppress_warnings==0)
990 >   fprintf(stderr, "Warning: skyclearness or skybrightness out of range in function direct_n_effi_PEREZ \n");*/
991  
992  
993   /* initialize category bounds (clearness index bounds) */
# Line 1062 | Line 1061 | return(value);
1061   /*check the range of epsilon and delta indexes of the perez parametrization*/
1062   void check_parametrization()
1063   {
1064 +
1065   if (skyclearness<skyclearinf || skyclearness>skyclearsup || skybrightness<skybriginf || skybrightness>skybrigsup)
1066                  {
1067  
1068   /*  limit sky clearness or sky brightness, 2009 11 13 by J. Wienold */
1069 +                
1070                  if (skyclearness<skyclearinf){
1071 <                        if (suppress_warnings==0)
1072 <                            /* fprintf(stderr,"Range warning: sky clearness too low (%lf)\n", skyclearness); */
1071 >                        /* if (suppress_warnings==0)
1072 >                            fprintf(stderr,"Range warning: sky clearness too low (%lf)\n", skyclearness); */
1073                          skyclearness=skyclearinf;
1074                  }
1075                  if (skyclearness>skyclearsup){
1076 <                        if (suppress_warnings==0)
1077 <                            /* fprintf(stderr,"Range warning: sky clearness too high (%lf)\n", skyclearness); */
1078 <                        skyclearness=skyclearsup-0.1;
1076 >                        /* if (suppress_warnings==0)
1077 >                            fprintf(stderr,"Range warning: sky clearness too high (%lf)\n", skyclearness); */
1078 >                        skyclearness=skyclearsup-0.001;
1079                  }
1080                  if (skybrightness<skybriginf){
1081 <                        if (suppress_warnings==0)
1082 <                            /* fprintf(stderr,"Range warning: sky brightness too low (%lf)\n", skybrightness); */
1081 >                        /* if (suppress_warnings==0)
1082 >                            fprintf(stderr,"Range warning: sky brightness too low (%lf)\n", skybrightness); */
1083                          skybrightness=skybriginf;
1084                  }
1085                  if (skybrightness>skybrigsup){
1086 <                        if (suppress_warnings==0)
1087 <                            /* fprintf(stderr,"Range warning: sky brightness too high (%lf)\n", skybrightness); */
1086 >                        /* if (suppress_warnings==0)
1087 >                            fprintf(stderr,"Range warning: sky brightness too high (%lf)\n", skybrightness); */
1088                          skybrightness=skybrigsup;
1089                  }
1090  
# Line 1092 | Line 1093 | if (skyclearness<skyclearinf || skyclearness>skyclears
1093   }
1094  
1095  
1096 +
1097 +
1098 +
1099   /* validity of the direct and diffuse components */
1100   void    check_illuminances()
1101   {
1102          if (directilluminance < 0) {
1103 <                fprintf(stderr,"WARNING: direct illuminance < 0. Using 0.0\n");
1103 >                if(suppress_warnings==0)
1104 >                { fprintf(stderr,"Warning: direct illuminance < 0. Using 0.0\n"); }
1105                  directilluminance = 0.0;
1106          }
1107          if (diffuseilluminance < 0) {
1108 <                fprintf(stderr,"WARNING: diffuse illuminance < 0. Using 0.0\n");
1108 >                if(suppress_warnings==0)
1109 >                { fprintf(stderr,"Warning: diffuse illuminance < 0. Using 0.0\n"); }
1110                  diffuseilluminance = 0.0;
1111          }
1112 <        if (directilluminance > solar_constant_l*1000.0) {
1113 <                fprintf(stderr,"ERROR: direct illuminance exceeds solar constant\n");
1114 <                exit(1);
1112 >        
1113 >        if (directilluminance+diffuseilluminance==0 && altitude > 0) {
1114 >                if(suppress_warnings==0)
1115 >                { fprintf(stderr,"Warning: zero illuminance at sun altitude > 0, printing error sky\n"); }
1116 >                print_error_sky();
1117 >                exit(0);
1118          }
1119 +        
1120 +        if (directilluminance > solar_constant_l) {
1121 +                if(suppress_warnings==0)
1122 +                { fprintf(stderr,"Warning: direct illuminance exceeds solar constant\n"); }
1123 +                print_error_sky();
1124 +                exit(0);
1125 +        }
1126   }
1127  
1128  
1129   void    check_irradiances()
1130   {
1131          if (directirradiance < 0) {
1132 <                fprintf(stderr,"WARNING: direct irradiance < 0. Using 0.0\n");
1132 >                if(suppress_warnings==0)
1133 >                { fprintf(stderr,"Warning: direct irradiance < 0. Using 0.0\n"); }
1134                  directirradiance = 0.0;
1135          }
1136          if (diffuseirradiance < 0) {
1137 <                fprintf(stderr,"WARNING: diffuse irradiance < 0. Using 0.0\n");
1137 >                if(suppress_warnings==0)
1138 >                { fprintf(stderr,"Warning: diffuse irradiance < 0. Using 0.0\n"); }
1139                  diffuseirradiance = 0.0;
1140          }
1141 +        
1142 +        if (directirradiance+diffuseirradiance==0 && altitude > 0) {
1143 +                if(suppress_warnings==0)
1144 +                { fprintf(stderr,"Warning: zero irradiance at sun altitude > 0, printing error sky\n"); }
1145 +                print_error_sky();
1146 +                exit(0);
1147 +        }
1148 +        
1149          if (directirradiance > solar_constant_e) {
1150 <                fprintf(stderr,"ERROR: direct irradiance exceeds solar constant\n");
1151 <                exit(1);
1150 >                if(suppress_warnings==0)
1151 >                { fprintf(stderr,"Warning: direct irradiance exceeds solar constant\n"); }
1152 >                print_error_sky();
1153 >                exit(0);
1154          }
1155   }
1156          
# Line 1181 | Line 1209 | void illu_to_irra_index()
1209   double  test1=0.1, test2=0.1, d_eff;
1210   int     counter=0;      
1211  
1212 < diffuseirradiance = diffuseilluminance*solar_constant_e/(solar_constant_l*1000);
1213 < directirradiance = directilluminance*solar_constant_e/(solar_constant_l*1000);
1212 > diffuseirradiance = diffuseilluminance*solar_constant_e/(solar_constant_l);
1213 > directirradiance = directilluminance*solar_constant_e/(solar_constant_l);
1214   skyclearness =  sky_clearness();
1215   skybrightness = sky_brightness();
1216   check_parametrization();
# Line 1210 | Line 1238 | while ( ((fabs(diffuseirradiance-test1)>10) || (fabs(d
1238          skybrightness = sky_brightness();
1239          skyclearness =  sky_clearness();
1240          check_parametrization();
1213 –        
1214 –        /*fprintf(stderr,"skyclearness = %lf, skybrightness = %lf, directirradiance = %lf, diffuseirradiance = %lf\n",skyclearness, skybrightness, directirradiance, diffuseirradiance);*/
1241                  
1242          }
1243  
# Line 1248 | Line 1274 | double calc_rel_lum_perez(double dzeta,double gamma,do
1274  
1275          if ( (epsilon <  skyclearinf) || (epsilon >= skyclearsup) )
1276          {
1277 <                fprintf(stderr,"Epsilon out of range in function calc_rel_lum_perez!\n");
1277 >                fprintf(stderr,"Error: epsilon out of range in function calc_rel_lum_perez!\n");
1278                  exit(1);
1279          }
1280  
# Line 1299 | Line 1325 | void coeff_lum_perez(double Z, double epsilon, double
1325  
1326          if ( (epsilon <  skyclearinf) || (epsilon >= skyclearsup) )
1327          {
1328 <                fprintf(stderr,"Epsilon out of range in function calc_rel_lum_perez !\n");
1328 >                fprintf(stderr,"Error: epsilon out of range in function coeff_lum_perez!\n");
1329                  exit(1);
1330          }
1331  
# Line 1368 | Line 1394 | void theta_phi_to_dzeta_gamma(double theta,double phi,
1394          else if ( (cos(Z)*cos(theta)+sin(Z)*sin(theta)*cos(phi)) > 1.1 )
1395          {
1396                  printf("error in calculation of gamma (angle between point and sun");
1397 <                exit(3);
1397 >                exit(1);
1398          }
1399          else
1400                  *gamma = acos(cos(Z)*cos(theta)+sin(Z)*sin(theta)*cos(phi));
# Line 1412 | Line 1438 | double         air_mass()
1438   double  m;
1439   if (sunzenith>90)
1440          {
1441 <        fprintf(stderr, "Solar zenith angle larger than 90 degrees in function air_mass()\n");
1442 <        exit(1);
1441 >        if(suppress_warnings==0)
1442 >        { fprintf(stderr, "Warning: air mass has reached the maximal value\n"); }
1443 >        sunzenith=90;
1444          }
1445   m = 1/( cos(sunzenith*M_PI/180)+0.15*exp( log(93.885-sunzenith)*(-1.253) ) );
1446   return(m);

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines