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

Comparing ray/src/gen/gendaylit.c (file contents):
Revision 2.13 by greg, Wed Aug 14 17:11:43 2013 UTC vs.
Revision 2.16 by greg, Wed Jul 30 17:30:27 2014 UTC

# Line 6 | Line 6
6   *                              1 Avenue Eugene Freyssinet, Saint-Quentin-Yvelines, France
7   */
8  
9 < #define _USE_MATH_DEFINES
9 > #define  _USE_MATH_DEFINES
10   #include  <stdio.h>
11   #include  <string.h>
12   #include  <math.h>
# Line 56 | 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 89 | 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 98 | 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
107
108 /* sun calculation constants */
109 extern double   s_latitude;
110 extern double   s_longitude;
111 extern double   s_meridian;
112
106   const double    AU = 149597890E3;
107   const double    solar_constant_e = 1367;    /* solar constant W/m^2 */
108 < const double    solar_constant_l = 127.5;   /* solar constant klux */
108 > const double    solar_constant_l = 127500;   /* solar constant lux */
109  
110   const double    half_sun_angle = 0.2665;
111   const double    half_direct_angle = 2.85;
112  
113 < const double    skyclearinf = 1.0;      /* limitations for the variation of the Perez parameters */
114 < const double    skyclearsup = 12.1;
113 > const double    skyclearinf = 1.0;          /* limitations for the variation of the Perez parameters */
114 > const double    skyclearsup = 12.01;
115   const double    skybriginf = 0.01;
116   const double    skybrigsup = 0.6;
117  
# Line 169 | Line 162 | float  timeinterval = 0;
162   char    *progname;
163   char    errmsg[128];
164  
165 + double  st;
166  
167  
174
168   int main(int argc, char** argv)
169   {
170          int  i;
# Line 182 | Line 175 | int main(int argc, char** argv)
175                  return 0;
176          }
177          if (argc < 4)
178 <                userror("arg count");
178 >                usage_error("arg count");
179          if (!strcmp(argv[1], "-ang")) {
180                  altitude = atof(argv[2]) * (M_PI/180);
181                  azimuth = atof(argv[3]) * (M_PI/180);
# Line 190 | Line 183 | int main(int argc, char** argv)
183          } else {
184                  month = atoi(argv[1]);
185                  if (month < 1 || month > 12)
186 <                        userror("bad month");
186 >                        usage_error("bad month");
187                  day = atoi(argv[2]);
188                  if (day < 1 || day > 31)
189 <                        userror("bad day");
189 >                        usage_error("bad day");
190                  hour = atof(argv[3]);
191                  if (hour < 0 || hour >= 24)
192 <                        userror("bad hour");
192 >                        usage_error("bad hour");
193                  tsolar = argv[3][0] == '+';
194          }
195          for (i = 4; i < argc; i++)
# Line 237 | Line 230 | int main(int argc, char** argv)
230                                  break;
231                          
232                          case 'O':
233 <                                output = atof(argv[++i]);       /*define the unit of the output of the program :
234 <                                                                sky and sun luminance/radiance (0==W visible, 1==W solar radiation, 2==lm) */
233 >                                output = atof(argv[++i]);       /*define the unit of the output of the program:
234 >                                                                sky and sun luminance/radiance
235 >                                                                (0==W visible, 1==W solar radiation, 2==lm) */
236                                  break;
237                                  
238                          case 'P':
# Line 270 | Line 264 | int main(int argc, char** argv)
264                                  globalirradiance = atof(argv[++i]);
265                                  break;
266                          
273                        /*
274                        case 'l':
275                                sunaltitude_border = atof(argv[++i]);
276                                break;
277                        */
278                                
267                          case 'i':
268                                  timeinterval = atof(argv[++i]);
269                                  break;
# Line 283 | Line 271 | int main(int argc, char** argv)
271                          
272                          default:
273                                  sprintf(errmsg, "unknown option: %s", argv[i]);
274 <                                userror(errmsg);
274 >                                usage_error(errmsg);
275                          }
276                  else
277 <                        userror("bad option");
277 >                        usage_error("bad option");
278  
279 <        if (fabs(s_meridian-s_longitude) > 30*M_PI/180)
280 <                fprintf(stderr,
293 <                    "%s: warning: %.1f hours btwn. standard meridian and longitude\n",
279 >        if (month && !tsolar && fabs(s_meridian-s_longitude) > 45*M_PI/180)
280 >                fprintf(stderr,"%s: warning: %.1f hours btwn. standard meridian and longitude\n",
281                      progname, (s_longitude-s_meridian)*12/M_PI);
282  
283  
284          /* dynamic memory allocation for the pointers */
298        
285          if ( (c_perez = calloc(5, sizeof(double))) == NULL )
286 <        {
301 <                fprintf(stderr,"Out of memory error in function main");
302 <                return 1;
303 <        }
286 >        { fprintf(stderr,"Out of memory error in function main"); return 1; }
287  
288 +        
289          printhead(argc, argv);
290          computesky();
307        
308        if(*(c_perez+1)>0)
309        {
310          fprintf(stderr, "Warning: positive Perez parameter B (= %lf), printing error sky\n",*(c_perez+1));      
311          print_error_sky();
312          exit(1);
313        }
314        
291          printsky();
292          return 0;
293 +
294   }
295  
296  
297  
298  
299  
323
324
325
300   void computesky()
301   {
302  
# Line 339 | Line 313 | void computesky()
313                  
314          if (month) {                    /* from date and time */
315                  int  jd;
316 <                double  sd, st;
316 >                double  sd;
317  
318                  jd = jdate(month, day);         /* Julian date */
319                  sd = sdec(jd);                  /* solar declination */
# Line 347 | Line 321 | void computesky()
321                          st = hour;
322                  else
323                          st = hour + stadj(jd);
350                                
324                  
352                if(st<solar_sunrise(month,day) || st>solar_sunset(month,day)) {
353                print_error_sky();
354                exit(1);
355                }
356                
325                                          
326                  if(timeinterval) {
327                          
# Line 361 | Line 329 | void computesky()
329                          fprintf(stderr, "time interval negative\n");
330                          exit(1);
331                          }
332 <                                        
333 <                        if(fabs(solar_sunrise(month,day)-st)<timeinterval/60) {                
334 <                        
335 <                        fprintf(stderr, "Solar position corrected at %d %d %.3f\n",month,day,hour);
336 <                        st= (st+timeinterval/120+solar_sunrise(month,day))/2;
332 >                                                                        
333 >                        if(fabs(solar_sunrise(month,day)-st)<=timeinterval/120) {                      
334 >                         st= (st+timeinterval/120+solar_sunrise(month,day))/2;
335 >                         if(suppress_warnings==0)
336 >                         { fprintf(stderr, "Solar position corrected at time step %d %d %.3f\n",month,day,hour); }
337                          }
338                  
339 <                        if(fabs(solar_sunset(month,day)-st)<timeinterval/60) {
340 <                        fprintf(stderr, "Solar position corrected at %d %d %.3f\n",month,day,hour);
341 <                        st= (st-timeinterval/120+solar_sunset(month,day))/2;
339 >                        if(fabs(solar_sunset(month,day)-st)<timeinterval/120) {
340 >                         st= (st-timeinterval/120+solar_sunset(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((st<solar_sunrise(month,day)-timeinterval/120) || (st>solar_sunset(month,day)+timeinterval/120)) {
346 +                          if(suppress_warnings==0)
347 +                          { fprintf(stderr, "Warning: sun position too low, printing error sky at %d %d %.3f\n",month,day,hour); }
348 +                         altitude = salt(sd, st);
349 +                         azimuth = sazi(sd, st);
350 +                         print_error_sky();
351 +                         exit(0);
352 +                        }
353                  }
354 +                else
355                  
356 +                if(st<solar_sunrise(month,day) || st>solar_sunset(month,day)) {
357 +                        if(suppress_warnings==0)
358 +                        { fprintf(stderr, "Warning: sun altitude below zero at time step %i %i %.2f, printing error sky\n",month,day,hour); }
359 +                        altitude = salt(sd, st);
360 +                        azimuth = sazi(sd, st);
361 +                        print_error_sky();
362 +                        exit(0);
363 +                }
364                  
365                  altitude = salt(sd, st);
366                  azimuth = sazi(sd, st);
# Line 385 | Line 372 | void computesky()
372          
373  
374  
388        /* if loop for the -l option. W.Sprenger (01/2013) */
389        /*
390        if (altitude*180/M_PI < sunaltitude_border) {
391        
392        if (suppress_warnings==0) {
393            fprintf(stderr, "Warning: sun altitude (%.3f degrees) below the border (%.3f degrees)\n",altitude*180/M_PI,sunaltitude_border);
394                                   }
395        print_error_sky();
396        exit(1);
397        }
398        */
399
375                          
376          if (!cloudy && altitude > 87.*M_PI/180.) {
377                  
# Line 408 | Line 383 | void computesky()
383                  altitude = 87.*M_PI/180.;
384          }
385          
386 +        
387 +        
388          sundir[0] = -sin(azimuth)*cos(altitude);
389          sundir[1] = -cos(azimuth)*cos(altitude);
390          sundir[2] = sin(altitude);
# Line 415 | Line 392 | void computesky()
392                  
393          /* calculation for the new functions */
394          sunzenith = 90 - altitude*180/M_PI;
395 <        
419 <        
395 >                        
396  
397          /* compute the inputs for the calculation of the light distribution over the sky*/
398          if (input==0)           /* P */
# Line 440 | Line 416 | void computesky()
416                  check_irradiances();
417                  skybrightness = sky_brightness();
418                  skyclearness =  sky_clearness();
419 +                
420                  check_parametrization();
421 <
421 >                                                        
422                  if (output==0 || output==2)
423                          {
424                          diffuseilluminance = diffuseirradiance*glob_h_diffuse_effi_PEREZ();/*diffuse horizontal illuminance*/
# Line 465 | Line 442 | void computesky()
442                          if (altitude<=0)
443                          {
444                                  if (suppress_warnings==0)
445 <                                     fprintf(stderr, "Warning: solar zenith angle larger than 90 degrees; using zero irradiance to proceed\n");
445 >                                     fprintf(stderr, "Warning: sun altitude < 0, proceed with irradiance values of zero\n");
446                                  directirradiance = 0;
447                                  diffuseirradiance = 0;
448                          } else {
# Line 505 | Line 482 | void computesky()
482                          {
483                                  if (suppress_warnings==0)
484                                          fprintf(stderr, "Warning: global irradiance is higher than the time-dependent solar constant s0\n");
508                                
485                                  globalirradiance=erbs_s0*0.999;                
486                          }
487                          
# Line 518 | Line 494 | void computesky()
494                          directirradiance=globalirradiance-diffuseirradiance;
495                          
496                          printf("# erbs_s0, erbs_kt, irr_dir_h, irr_diff: %.3f %.3f %.3f %.3f\n", erbs_s0, erbs_kt, directirradiance, diffuseirradiance);
497 <                        printf("# WARNING: the -E option is only recommended for a rough estimation!");
497 >                        printf("# WARNING: the -E option is only recommended for a rough estimation!\n");
498                          
499                          directirradiance=directirradiance/sin(altitude);
500                                                                                                                          
# Line 541 | Line 517 | void computesky()
517                  
518                  
519          
520 <        else    {fprintf(stderr,"error at the input arguments"); exit(1);}
520 >        else    { fprintf(stderr,"error at the input arguments"); exit(1); }
521  
522  
523          
# Line 563 | Line 539 | void computesky()
539  
540          
541          
566        
567        
542          /*calculation of the modelled luminance */
543          for (j=0;j<145;j++)
544          {
# Line 610 | Line 584 | void computesky()
584          else
585                  solarradiance = directilluminance/(2*M_PI*(1-cos(half_sun_angle*M_PI/180)));
586          
613                        
587  
588  
589 < /* Compute the ground radiance */
590 < zenithbr=calc_rel_lum_perez(0.0,radians(sunzenith),radians(sunzenith),skyclearness,skybrightness,coeff_perez);
591 < zenithbr*=diffnormalization;
589 >        /* Compute the ground radiance */
590 >        zenithbr=calc_rel_lum_perez(0.0,radians(sunzenith),radians(sunzenith),skyclearness,skybrightness,coeff_perez);
591 >        zenithbr*=diffnormalization;
592          
593 < if (skyclearness==1)
593 >        if (skyclearness==1)
594          normfactor = 0.777778;
595                  
596 < if (skyclearness>=6)
596 >        if (skyclearness>=6)
597          {              
598          F2 = 0.274*(0.91 + 10.0*exp(-3.0*(M_PI/2.0-altitude)) + 0.45*sundir[2]*sundir[2]);
599          normfactor = normsc()/F2/M_PI;
600          }
601  
602 < if ( (skyclearness>1) && (skyclearness<6) )
602 >        if ( (skyclearness>1) && (skyclearness<6) )
603          {
604          S_INTER=1;
605          F2 = (2.739 + .9891*sin(.3119+2.6*altitude)) * exp(-(M_PI/2.0-altitude)*(.4441+1.48*altitude));
606          normfactor = normsc()/F2/M_PI;
607          }
608  
609 < groundbr = zenithbr*normfactor;
609 >        groundbr = zenithbr*normfactor;
610  
611 < if (dosun&&(skyclearness>1))
611 >        if (dosun&&(skyclearness>1))
612          groundbr += 6.8e-5/M_PI*solarradiance*sundir[2];                
613  
614 < groundbr *= gprefl;
614 >        groundbr *= gprefl;
615  
616  
617 +                
618 +        if(*(c_perez+1)>0)
619 +        {
620 +          if(suppress_warnings==0)
621 +                {  fprintf(stderr, "Warning: positive Perez parameter B (= %lf), printing error sky\n",*(c_perez+1));}  
622 +          print_error_sky();
623 +          exit(0);
624 +        }
625  
626 +
627   return;
628   }
629  
630  
631  
632  
651 void print_error_sky()
652 {
653        sundir[0] = -sin(azimuth)*cos(altitude);
654        sundir[1] = -cos(azimuth)*cos(altitude);
655        sundir[2] = sin(altitude);
656        
657        printf("\nvoid brightfunc skyfunc\n");
658        printf("2 skybright perezlum.cal\n");
659        printf("0\n");
660        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]);
661 }
662        
633  
664
665
666
667
668
634   double solar_sunset(int month,int day)
635   {
636       float W;
# Line 675 | Line 640 | double solar_sunset(int month,int day)
640   }
641  
642  
643 +
644 +
645   double solar_sunrise(int month,int day)
646   {
647       float W;
# Line 686 | Line 653 | double solar_sunrise(int month,int day)
653  
654  
655  
656 + void printsky()
657 + {      
658 +        
659 +        printf("# Local solar time: %.2f\n", st);
660 +        printf("# Solar altitude and azimuth: %.1f %.1f\n", altitude*180/M_PI, azimuth*180/M_PI);
661  
662  
691
692
693
694
695
696 void printsky()                 /* print out sky */
697 {
663          if (dosun&&(skyclearness>1))
664          {              
665                  printf("\nvoid light solar\n");
# Line 723 | Line 688 | void printsky()                        /* print out sky */
688   }
689  
690  
691 +
692 + void print_error_sky()
693 + {
694 +
695 +
696 +        sundir[0] = -sin(azimuth)*cos(altitude);
697 +        sundir[1] = -cos(azimuth)*cos(altitude);
698 +        sundir[2] = sin(altitude);
699 +
700 +        printf("# Local solar time: %.2f\n", st);
701 +        printf("# Solar altitude and azimuth: %.1f %.1f\n", altitude*180/M_PI, azimuth*180/M_PI);
702 +
703 +        printf("\nvoid brightfunc skyfunc\n");
704 +        printf("2 skybright perezlum.cal\n");
705 +        printf("0\n");
706 +        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]);
707 + }
708 +        
709 +
710 +
711 +
712 +
713   void printdefaults()                    /* print default values */
714   {
715          printf("-g %f\t\t\t# Ground plane reflectance\n", gprefl);
# Line 736 | Line 723 | void printdefaults()                   /* print default values */
723   }
724  
725  
726 < void userror(char* msg)                 /* print usage error and quit */
726 >
727 >
728 > void usage_error(char* msg)                     /* print usage error and quit */
729   {
730          if (msg != NULL)
731                  fprintf(stderr, "%s: Use error - %s\n\n", progname, msg);
# Line 753 | Line 742 | void userror(char* msg)                        /* print usage error and qui
742          fprintf(stderr, "       -E global-horizontal-irradiance (W/m^2)\n\n");
743          fprintf(stderr, "       Output specification with option:\n");
744          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");
745 <        fprintf(stderr, "       gendaylit version 2.3 (2013/08/08)  \n\n");
745 >        fprintf(stderr, "       gendaylit version 2.4 (2013/09/04)  \n\n");
746          exit(1);
747   }
748  
749  
750  
751 +
752   double normsc()           /* compute normalization factor (E0*F2/L0) */
753   {
754          static double  nfc[2][5] = {
# Line 782 | Line 772 | double normsc()                  /* compute normalization fac
772  
773  
774  
775 +
776 +
777   void printhead(int ac, char** av)               /* print command header */
778   {
779          putchar('#');
# Line 795 | Line 787 | void printhead(int ac, char** av)              /* print command he
787  
788  
789  
790 +
791 +
792   /* Perez models */
793  
794   /* Perez global horizontal luminous efficacy model */
# Line 866 | Line 860 | double glob_h_effi_PEREZ()
860  
861  
862  
869
863          for (i=1; i<=category_total_number; i++)
864          {
865                  if ( (skyclearness >= category_bounds[i]) && (skyclearness < category_bounds[i+1]) )
# Line 880 | Line 873 | double glob_h_effi_PEREZ()
873   }
874  
875  
876 +
877 +
878   /* global horizontal diffuse efficacy model, according to PEREZ */
879   double glob_h_diffuse_effi_PEREZ()
880   {
# Line 887 | Line 882 | double glob_h_diffuse_effi_PEREZ()
882          double  category_bounds[10], a[10], b[10], c[10], d[10];
883          int     category_total_number, category_number, i;
884  
890
891
892        
885          check_parametrization();
886          
887          
888   /*if ((skyclearness<skyclearinf || skyclearness>skyclearsup || skybrightness<skybriginf || skybrightness>skybrigsup) && suppress_warnings==0)
889 <     fprintf(stderr, "Warning: skyclearness or skybrightness out of range in function glob_h_diffuse_PEREZ \n"); */
889 > fprintf(stderr, "Warning: skyclearness or skybrightness out of range in function glob_h_diffuse_PEREZ \n"); */
890      
899
900
891   /* initialize category bounds (clearness index bounds) */
892  
893          category_total_number = 8;
# Line 953 | Line 943 | double glob_h_diffuse_effi_PEREZ()
943  
944  
945  
956
946          category_number = -1;
947          for (i=1; i<=category_total_number; i++)
948          {
# Line 963 | Line 952 | double glob_h_diffuse_effi_PEREZ()
952  
953          if (category_number == -1) {
954                  if (suppress_warnings==0)
955 <                    fprintf(stderr, "ERROR: Model parameters out of range, skyclearness = %lf \n", skyclearness);
955 >                    fprintf(stderr, "Warning: sky clearness (= %.3f) too high, printing error sky\n", skyclearness);
956                  print_error_sky();
957 <                exit(1);
957 >                exit(0);
958          }
959                  
960  
# Line 978 | Line 967 | double glob_h_diffuse_effi_PEREZ()
967  
968  
969  
970 +
971 +
972 +
973   /* direct normal efficacy model, according to PEREZ */
974  
975   double direct_n_effi_PEREZ()
# Line 988 | Line 980 | double         category_bounds[10], a[10], b[10], c[10], d[10
980   int     category_total_number, category_number, i;
981  
982  
983 < if ((skyclearness<skyclearinf || skyclearness>skyclearsup || skybrightness<skybriginf || skybrightness>skybrigsup) && suppress_warnings==0)
984 <   fprintf(stderr, "Warning: skyclearness or skybrightness out of range in function direct_n_effi_PEREZ \n");
983 > /*if ((skyclearness<skyclearinf || skyclearness>skyclearsup || skybrightness<skybriginf || skybrightness>skybrigsup) && suppress_warnings==0)
984 >   fprintf(stderr, "Warning: skyclearness or skybrightness out of range in function direct_n_effi_PEREZ \n");*/
985  
986  
987   /* initialize category bounds (clearness index bounds) */
# Line 1063 | Line 1055 | return(value);
1055   /*check the range of epsilon and delta indexes of the perez parametrization*/
1056   void check_parametrization()
1057   {
1058 +
1059   if (skyclearness<skyclearinf || skyclearness>skyclearsup || skybrightness<skybriginf || skybrightness>skybrigsup)
1060                  {
1061  
1062   /*  limit sky clearness or sky brightness, 2009 11 13 by J. Wienold */
1063 +                
1064                  if (skyclearness<skyclearinf){
1065 <                        if (suppress_warnings==0)
1066 <                            /* fprintf(stderr,"Range warning: sky clearness too low (%lf)\n", skyclearness); */
1065 >                        /* if (suppress_warnings==0)
1066 >                            fprintf(stderr,"Range warning: sky clearness too low (%lf)\n", skyclearness); */
1067                          skyclearness=skyclearinf;
1068                  }
1069                  if (skyclearness>skyclearsup){
1070 <                        if (suppress_warnings==0)
1071 <                            /* fprintf(stderr,"Range warning: sky clearness too high (%lf)\n", skyclearness); */
1072 <                        skyclearness=skyclearsup-0.1;
1070 >                        /* if (suppress_warnings==0)
1071 >                            fprintf(stderr,"Range warning: sky clearness too high (%lf)\n", skyclearness); */
1072 >                        skyclearness=skyclearsup-0.001;
1073                  }
1074                  if (skybrightness<skybriginf){
1075 <                        if (suppress_warnings==0)
1076 <                            /* fprintf(stderr,"Range warning: sky brightness too low (%lf)\n", skybrightness); */
1075 >                        /* if (suppress_warnings==0)
1076 >                            fprintf(stderr,"Range warning: sky brightness too low (%lf)\n", skybrightness); */
1077                          skybrightness=skybriginf;
1078                  }
1079                  if (skybrightness>skybrigsup){
1080 <                        if (suppress_warnings==0)
1081 <                            /* fprintf(stderr,"Range warning: sky brightness too high (%lf)\n", skybrightness); */
1080 >                        /* if (suppress_warnings==0)
1081 >                            fprintf(stderr,"Range warning: sky brightness too high (%lf)\n", skybrightness); */
1082                          skybrightness=skybrigsup;
1083                  }
1084  
# Line 1093 | Line 1087 | if (skyclearness<skyclearinf || skyclearness>skyclears
1087   }
1088  
1089  
1090 +
1091 +
1092 +
1093   /* validity of the direct and diffuse components */
1094   void    check_illuminances()
1095   {
1096          if (directilluminance < 0) {
1097 <                fprintf(stderr,"WARNING: direct illuminance < 0. Using 0.0\n");
1097 >                if(suppress_warnings==0)
1098 >                { fprintf(stderr,"Warning: direct illuminance < 0. Using 0.0\n"); }
1099                  directilluminance = 0.0;
1100          }
1101          if (diffuseilluminance < 0) {
1102 <                fprintf(stderr,"WARNING: diffuse illuminance < 0. Using 0.0\n");
1102 >                if(suppress_warnings==0)
1103 >                { fprintf(stderr,"Warning: diffuse illuminance < 0. Using 0.0\n"); }
1104                  diffuseilluminance = 0.0;
1105          }
1106 <        if (directilluminance > solar_constant_l*1000.0) {
1107 <                fprintf(stderr,"ERROR: direct illuminance exceeds solar constant\n");
1108 <                exit(1);
1106 >        
1107 >        if (directilluminance+diffuseilluminance==0 && altitude > 0) {
1108 >                if(suppress_warnings==0)
1109 >                { fprintf(stderr,"Warning: zero illuminance at sun altitude > 0, printing error sky\n"); }
1110 >                print_error_sky();
1111 >                exit(0);
1112          }
1113 +        
1114 +        if (directilluminance > solar_constant_l) {
1115 +                if(suppress_warnings==0)
1116 +                { fprintf(stderr,"Warning: direct illuminance exceeds solar constant\n"); }
1117 +                print_error_sky();
1118 +                exit(0);
1119 +        }
1120   }
1121  
1122  
1123   void    check_irradiances()
1124   {
1125          if (directirradiance < 0) {
1126 <                fprintf(stderr,"WARNING: direct irradiance < 0. Using 0.0\n");
1126 >                if(suppress_warnings==0)
1127 >                { fprintf(stderr,"Warning: direct irradiance < 0. Using 0.0\n"); }
1128                  directirradiance = 0.0;
1129          }
1130          if (diffuseirradiance < 0) {
1131 <                fprintf(stderr,"WARNING: diffuse irradiance < 0. Using 0.0\n");
1131 >                if(suppress_warnings==0)
1132 >                { fprintf(stderr,"Warning: diffuse irradiance < 0. Using 0.0\n"); }
1133                  diffuseirradiance = 0.0;
1134          }
1135 +        
1136 +        if (directirradiance+diffuseirradiance==0 && altitude > 0) {
1137 +                if(suppress_warnings==0)
1138 +                { fprintf(stderr,"Warning: zero irradiance at sun altitude > 0, printing error sky\n"); }
1139 +                print_error_sky();
1140 +                exit(0);
1141 +        }
1142 +        
1143          if (directirradiance > solar_constant_e) {
1144 <                fprintf(stderr,"ERROR: direct irradiance exceeds solar constant\n");
1145 <                exit(1);
1144 >                if(suppress_warnings==0)
1145 >                { fprintf(stderr,"Warning: direct irradiance exceeds solar constant\n"); }
1146 >                print_error_sky();
1147 >                exit(0);
1148          }
1149   }
1150          
# Line 1182 | Line 1203 | void illu_to_irra_index()
1203   double  test1=0.1, test2=0.1, d_eff;
1204   int     counter=0;      
1205  
1206 < diffuseirradiance = diffuseilluminance*solar_constant_e/(solar_constant_l*1000);
1207 < directirradiance = directilluminance*solar_constant_e/(solar_constant_l*1000);
1206 > diffuseirradiance = diffuseilluminance*solar_constant_e/(solar_constant_l);
1207 > directirradiance = directilluminance*solar_constant_e/(solar_constant_l);
1208   skyclearness =  sky_clearness();
1209   skybrightness = sky_brightness();
1210   check_parametrization();
# Line 1211 | Line 1232 | while ( ((fabs(diffuseirradiance-test1)>10) || (fabs(d
1232          skybrightness = sky_brightness();
1233          skyclearness =  sky_clearness();
1234          check_parametrization();
1214        
1215        /*fprintf(stderr,"skyclearness = %lf, skybrightness = %lf, directirradiance = %lf, diffuseirradiance = %lf\n",skyclearness, skybrightness, directirradiance, diffuseirradiance);*/
1235                  
1236          }
1237  
# Line 1249 | Line 1268 | double calc_rel_lum_perez(double dzeta,double gamma,do
1268  
1269          if ( (epsilon <  skyclearinf) || (epsilon >= skyclearsup) )
1270          {
1271 <                fprintf(stderr,"Epsilon out of range in function calc_rel_lum_perez!\n");
1271 >                fprintf(stderr,"Error: epsilon out of range in function calc_rel_lum_perez!\n");
1272                  exit(1);
1273          }
1274  
# Line 1300 | Line 1319 | void coeff_lum_perez(double Z, double epsilon, double
1319  
1320          if ( (epsilon <  skyclearinf) || (epsilon >= skyclearsup) )
1321          {
1322 <                fprintf(stderr,"Epsilon out of range in function calc_rel_lum_perez !\n");
1322 >                fprintf(stderr,"Error: epsilon out of range in function coeff_lum_perez!\n");
1323                  exit(1);
1324          }
1325  
# Line 1369 | Line 1388 | void theta_phi_to_dzeta_gamma(double theta,double phi,
1388          else if ( (cos(Z)*cos(theta)+sin(Z)*sin(theta)*cos(phi)) > 1.1 )
1389          {
1390                  printf("error in calculation of gamma (angle between point and sun");
1391 <                exit(3);
1391 >                exit(1);
1392          }
1393          else
1394                  *gamma = acos(cos(Z)*cos(theta)+sin(Z)*sin(theta)*cos(phi));
# Line 1413 | Line 1432 | double         air_mass()
1432   double  m;
1433   if (sunzenith>90)
1434          {
1435 <        fprintf(stderr, "Solar zenith angle larger than 90 degrees in function air_mass()\n");
1436 <        exit(1);
1435 >        if(suppress_warnings==0)
1436 >        { fprintf(stderr, "Warning: air mass has reached the maximal value\n"); }
1437 >        sunzenith=90;
1438          }
1439   m = 1/( cos(sunzenith*M_PI/180)+0.15*exp( log(93.885-sunzenith)*(-1.253) ) );
1440   return(m);

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines