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.15 by greg, Mon Nov 18 18:07:16 2013 UTC vs.
Revision 2.21 by greg, Thu Jan 28 19:03:15 2021 UTC

# Line 4 | Line 4
4   *                              Centre de Valbonne, 500 route des Lucioles, 06565 Sophia Antipolis Cedex, France
5   *                              *BOUYGUES
6   *                              1 Avenue Eugene Freyssinet, Saint-Quentin-Yvelines, France
7 + *  print colored output if activated in command line (-C). Based on model from A. Diakite, TU-Berlin. Implemented by J. Wienold, August 26 2018
8 + *  version 2.6 (2021/01/29): dew point dependency added according to Perez publication 1990 (W -> atm_preci_water=exp(0.07*Td-0.075) ). by J. Wienold, EPFL
9   */
10  
11   #define  _USE_MATH_DEFINES
# Line 53 | Line 55 | float defangle_phi[] = {
55          90, 105, 120, 135, 150, 165, 180, 195, 210, 225, 240, 255, 270, 285, 300, 315, 330, 345, 0, 20, 40, 60,
56          80, 100, 120, 140, 160, 180, 200, 220, 240, 260, 280, 300, 320, 340, 0, 30, 60, 90, 120, 150, 180, 210,
57          240, 270, 300, 330, 0, 60, 120, 180, 240, 300, 0};
58 + /* default values for Berlin */
59 + float   locus[] = {
60 + -4.843e9,2.5568e6,0.24282e3,0.23258,-4.843e9,2.5568e6,0.24282e3,0.23258,-1.2848,1.7519,-0.093786};
61  
62  
63  
# Line 100 | Line 105 | double         air_mass();
105  
106   double  solar_sunset(int month, int day);
107   double  solar_sunrise(int month, int day);
103 – double  stadj();
104 – int     jdate(int month, int day);
108  
106 –
107 – /* sun calculation constants */
108 – extern double   s_latitude;
109 – extern double   s_longitude;
110 – extern double   s_meridian;
111 –
109   const double    AU = 149597890E3;
110   const double    solar_constant_e = 1367;    /* solar constant W/m^2 */
111   const double    solar_constant_l = 127500;   /* solar constant lux */
# Line 124 | Line 121 | const double   skybrigsup = 0.6;
121  
122  
123   /* required values */
124 + int  year = 0;                                  /* year (optional) */
125   int     month, day;                             /* date */
126   double  hour;                                   /* time */
127   int     tsolar;                                 /* 0=standard, 1=solar */
# Line 136 | Line 134 | double skyclearness = 0;
134   double  skybrightness = 0;
135   double  solarradiance;
136   double  diffuseilluminance, directilluminance, diffuseirradiance, directirradiance, globalirradiance;
137 < double  sunzenith, daynumber, atm_preci_water=2;
137 > double  sunzenith, daynumber, atm_preci_water, Td=10.97353115;
138  
139   /*double  sunaltitude_border = 0;*/
140   double  diffnormalization = 0;
# Line 146 | Line 144 | double         *c_perez;
144   int     output=0;       /* define the unit of the output (sky luminance or radiance): */
145                          /* visible watt=0, solar watt=1, lumen=2 */
146   int     input=0;        /* define the input for the calulation */
147 <
147 > int     color_output=0;
148   int     suppress_warnings=0;
149  
150          /* default values */
# Line 157 | Line 155 | double  betaturbidity = 0.1;
155   double  gprefl = 0.2;
156   int     S_INTER=0;
157  
158 +
159          /* computed values */
160   double  sundir[3];
161   double  groundbr = 0;
# Line 201 | Line 200 | int main(int argc, char** argv)
200          for (i = 4; i < argc; i++)
201                  if (argv[i][0] == '-' || argv[i][0] == '+')
202                          switch (argv[i][1]) {
203 +                        case 'd':
204 +                                Td = atof(argv[++i]);
205 +                                if (Td < -40 || Td > 40) {
206 +                                        Td=10.97353115; }
207 +                                break;
208                          case 's':
209                                  cloudy = 0;
210                                  dosun = argv[i][0] == '+';
211                                  break;
212 +                        case 'y':
213 +                                year = atoi(argv[++i]);
214 +                                break;
215                          case 'R':
216                                  u_solar = argv[i][1] == 'R' ? -1 : 1;
217                                  solarbr = atof(argv[++i]);
# Line 213 | Line 220 | int main(int argc, char** argv)
220                                  cloudy = argv[i][0] == '+' ? 2 : 1;
221                                  dosun = 0;
222                                  break;
223 +                        case 'C':
224 +                                if (argv[i][2] == 'I' && argv[i][3] == 'E' ) {
225 +                                locus[0] = -4.607e9;
226 +                                locus[1] = 2.9678e6;
227 +                                locus[2] = 0.09911e3;
228 +                                locus[3] = 0.244063;
229 +                                locus[4] = -2.0064e9;
230 +                                locus[5] = 1.9018e6;
231 +                                locus[6] = 0.24748e3;
232 +                                locus[7] = 0.23704;
233 +                                locus[8] = -3.0;
234 +                                locus[9] = 2.87;
235 +                                locus[10] = -0.275;
236 +                                 }else{ color_output = 1;
237 +                                 }
238 +                                break;
239 +                        case 'l':
240 +                                locus[0] = atof(argv[++i]);
241 +                                locus[1] = atof(argv[++i]);
242 +                                locus[2] = atof(argv[++i]);
243 +                                locus[3] = atof(argv[++i]);
244 +                                locus[4] = locus[0];
245 +                                locus[5] = locus[1];
246 +                                locus[6] = locus[2];
247 +                                locus[7] = locus[3];
248 +                                locus[8] = atof(argv[++i]);
249 +                                locus[9] = atof(argv[++i]);
250 +                                locus[10] = atof(argv[++i]);
251 +                                break;
252 +
253                          case 't':
254                                  betaturbidity = atof(argv[++i]);
255                                  break;
# Line 292 | Line 329 | int main(int argc, char** argv)
329          { fprintf(stderr,"Out of memory error in function main"); return 1; }
330  
331          
332 +        atm_preci_water=exp(0.07*Td-0.075);
333          printhead(argc, argv);
334          computesky();
335          printsky();
# Line 318 | Line 356 | void computesky()
356          /* compute solar direction */
357                  
358          if (month) {                    /* from date and time */
321 –                int  jd;
359                  double  sd;
360  
361 <                jd = jdate(month, day);         /* Julian date */
362 <                sd = sdec(jd);                  /* solar declination */
363 <                if (tsolar)                     /* solar time */
364 <                        st = hour;
365 <                else
366 <                        st = hour + stadj(jd);
367 <                
368 <                                        
361 >                st = hour;
362 >                if (year) {                     /* Michalsky algorithm? */
363 >                        double  mjd = mjdate(year, month, day, hour);
364 >                        if (tsolar)
365 >                                sd = msdec(mjd, NULL);
366 >                        else
367 >                                sd = msdec(mjd, &st);
368 >                } else {
369 >                        int  jd = jdate(month, day);    /* Julian date */
370 >                        sd = sdec(jd);                  /* solar declination */
371 >                        if (!tsolar)                    /* get solar time? */
372 >                                st = hour + stadj(jd);
373 >                }
374 >                                                        
375                  if(timeinterval) {
376                          
377                          if(timeinterval<0) {
# Line 664 | Line 707 | void printsky()
707          
708          printf("# Local solar time: %.2f\n", st);
709          printf("# Solar altitude and azimuth: %.1f %.1f\n", altitude*180/M_PI, azimuth*180/M_PI);
710 +        printf("# epsilon, delta, atmospheric precipitable water content : %.4f %.4f %.4f \n", skyclearness, skybrightness,atm_preci_water );
711  
712  
713          if (dosun&&(skyclearness>1))
# Line 682 | Line 726 | void printsky()
726                  printf("0\n0\n");
727                  printf("4 %f %f %f %f\n", sundir[0], sundir[1], sundir[2], 2*half_sun_angle);
728          }
729 <        
730 <
731 <        printf("\nvoid brightfunc skyfunc\n");
732 <        printf("2 skybright perezlum.cal\n");
729 > /* print colored output if activated in command line (-C). Based on model from A. Diakite, TU-Berlin. Implemented by J. Wienold, August 26 2018 */      
730 >        if  (color_output==1 && skyclearness < 4.5 && skyclearness >1.065 )  
731 >        {
732 >        fprintf(stderr, "       warning: sky clearness(epsilon)= %f \n",skyclearness);
733 >        fprintf(stderr, "       warning: intermediate sky!! \n");
734 >        fprintf(stderr, "       warning: color model for intermediate sky pending  \n");
735 >        fprintf(stderr, "       warning: no color output ! \n");
736 >        color_output=0;
737 >        }
738 >        if (color_output==1)
739 >        {
740 >        printf("\nvoid colorfunc skyfunc\n");
741 >        printf("4 skybright_r skybright_g skybright_b perezlum_c.cal\n");
742          printf("0\n");
743 <        printf("10 %.3e %.3e %lf %lf %lf %lf %lf %f %f %f \n", diffnormalization, groundbr,
743 >        printf("22 %.3e %.3e %lf %lf %lf %lf %lf %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f\n", diffnormalization, groundbr,
744                  *(c_perez+0),*(c_perez+1),*(c_perez+2),*(c_perez+3),*(c_perez+4),
745 <                sundir[0], sundir[1], sundir[2]);
745 >                sundir[0], sundir[1], sundir[2],skyclearness,locus[0],locus[1],locus[2],locus[3],locus[4],locus[5],locus[6],locus[7],locus[8],locus[9],locus[10]);
746 >        }else{
747 >        printf("\nvoid brightfunc skyfunc\n");
748 >        printf("2 skybright perezlum.cal\n");
749 >        printf("0\n");
750 >        printf("10 %.3e %.3e %lf %lf %lf %lf %lf %f %f %f \n", diffnormalization, groundbr,
751 >                *(c_perez+0),*(c_perez+1),*(c_perez+2),*(c_perez+3),*(c_perez+4),
752 >                sundir[0], sundir[1], sundir[2]);
753 >         }
754          
755   }
756  
# Line 735 | Line 796 | void usage_error(char* msg)                    /* print usage error and
796   {
797          if (msg != NULL)
798                  fprintf(stderr, "%s: Use error - %s\n\n", progname, msg);
799 <        fprintf(stderr, "Usage: %s      month day hour    [...]\n", progname);
800 <        fprintf(stderr, "   or: %s -ang altitude azimuth  [...]\n", progname);
799 >        fprintf(stderr, "Usage: %s      month day hour [-y year]        [...]\n", progname);
800 >        fprintf(stderr, "   or: %s -ang altitude azimuth                [...]\n", progname);
801          fprintf(stderr, "               followed by:      -P          epsilon delta [options]\n");
802          fprintf(stderr, "                        or:      [-W|-L|-G]  direct_value diffuse_value [options]\n");
803 <        fprintf(stderr, "                        or:      -E          global_irradiance [options]\n\n");
804 <        fprintf(stderr, "       Description:\n");
803 >        fprintf(stderr, "                        or:      -E          global_irradiance [options]\n\n");
804 >        fprintf(stderr, "    Description:\n");
805          fprintf(stderr, "       -P epsilon delta  (these are the Perez parameters) \n");
806          fprintf(stderr, "       -W direct-normal-irradiance diffuse-horizontal-irradiance (W/m^2)\n");
807          fprintf(stderr, "       -L direct-normal-illuminance diffuse-horizontal-illuminance (lux)\n");
# Line 748 | Line 809 | void usage_error(char* msg)                    /* print usage error and
809          fprintf(stderr, "       -E global-horizontal-irradiance (W/m^2)\n\n");
810          fprintf(stderr, "       Output specification with option:\n");
811          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");
812 <        fprintf(stderr, "       gendaylit version 2.4 (2013/09/04)  \n\n");
812 >        fprintf(stderr, "       gendaylit version 2.5 (2018/04/18)  \n\n");
813          exit(1);
814   }
815  
# Line 1374 | Line 1435 | void coeff_lum_perez(double Z, double epsilon, double
1435   /* degrees into radians */
1436   double radians(double degres)
1437   {
1438 <        return degres*M_PI/180.0;
1438 >        return degres*(M_PI/180.);
1439   }
1440  
1441  
1442   /* radian into degrees */
1443   double degres(double radians)
1444   {
1445 <        return radians/M_PI*180.0;
1445 >        return radians*(180./M_PI);
1446   }
1447  
1448  
# Line 1412 | Line 1473 | double integ_lv(float *lv,float *theta)
1473                  buffer += (*(lv+i))*cos(radians(*(theta+i)));
1474          }
1475                          
1476 <        return buffer*2*M_PI/144;
1476 >        return buffer*(2.*M_PI/145.);
1477   }
1478  
1479  

Diff Legend

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