/usr/share/radiance/Earth.cal is in radiance-materials 4R1+20120125-1.1.
This file is owned by root:root, with mode 0o644.
The actual contents of the file can be viewed below.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 | {
Coordinate mapping for point above Earth
assuming uniform projection of latitude and longitude
5/1/96 Greg Ward
A1 = latitude position (degrees)
A2 = longitude position (degrees)
A3 = height above Earth (kilometers)
A4 = Julian date (1-365)
A5 = Greenwich Mean Time (0-24)
By modifying this file so that the above arguments are defined
instead as constants, significant gains in speed are possible.
(Be sure to change Su, Sv and Sw to constant definitions also.)
Y-axis points North, Z-axis points up (away from Earth)
}
{
Example usage:
void colorpict epict
7 Eval Eval Eval Earth.pic Earth.cal u v
0
6 lat long height JDate GMT clouds
epict glow eglow
0
0
4 1 1 1 0
eglow source earth
0
0
4 0 0 -1 max_angle
where:
max_angle = 2*180/PI * asin(EarthRad/(EarthRad+height))
If you are really close to the Earth (below atmosphere), use:
!gensky month day time -a lat -o lon -m std_meridian
skybright colorpict epict
7 noop noop noop Earth.pic Earth.cal u v
0
3 lat long height
epict glow eglow
0
0
4 1 1 1 0
eglow source earth
0
0
4 0 0 -1 max_angle
}
EarthRad : 6378.5; { Earth radius in kilometers }
EarthTilt : 23.45*DEGREE; { Earth's orbital axis tilt }
EarthIrrad : 770 {watts/meter^2}; { Maximum Earth irradiance }
PA : 256/512; { Picture Height/Width }
CloudRho : 0.9; { Cloud reflectance }
{ Change the following to constants to speed things up }
Lat = arg(1)*DEGREE; { Position latitude }
Lon = arg(2)*DEGREE; { Position longitude }
Height = arg(3); { Position height above Earth (km) }
JDate = arg(4); { Julian date (1-365) }
GMT = arg(5); { Greenwich Mean Time (0-24) }
Cloudy = if(arg(0)-5.5, arg(6), 0);
GOu = EarthRad*sin(Lon)*cos(Lat); { Starting Earth position }
GOv = EarthRad*cos(Lon)*cos(Lat);
GOw = EarthRad*sin(Lat);
raduv = max(FTINY, sqrt(GOu*GOu + GOv*GOv)); { Earthwise direction vector }
GVXu = -GOv / raduv;
GVXv = GOu / raduv;
GVXw = 0;
GVYu = (-GOw*GOu) / EarthRad / raduv;
GVYv = (GOw*-GOv) / EarthRad / raduv;
GVYw = (GOu*GOu + GOv*GOv) / EarthRad / raduv;
radxy = max(FTINY, sqrt(Dx*Dx + Dy*Dy));
GDu = (GVXu*Dx + GVYu*Dy) / radxy;
GDv = (GVXv*Dx + GVYv*Dy) / radxy;
GDw = (GVXw*Dx + GVYw*Dy) / radxy;
tana = if(FTINY+Dz, radxy/FTINY, radxy/-Dz);
det = sq(EarthRad+Height) - (tana*tana + 1)*(2*EarthRad*Height + Height*Height);
a = if(tana-FTINY, b/tana, Height);
b = if(tana-FTINY, (EarthRad + Height - Sqrt(det)) / (tana + 1/tana), 0);
c = EarthRad + Height - a;
digfactor = 1 - (a - Height)/EarthRad;
{ Global position vector }
GPu = GOu*digfactor + GDu*b;
GPv = GOv*digfactor + GDv*b;
GPw = GOw*digfactor + GDw*b;
{ Earthly coordinates }
longitude = atan2(GPu, GPv);
latitude = Asin(GPw/EarthRad);
{ Image lookup coordinates }
u = if(PA-1, 1, 1/PA) * (.5 - longitude/(2*PI));
v = if(PA-1, PA, 1) * (.5 + latitude/PI);
{ Solar declination (lat.) }
Slat = EarthTilt * sin(2*PI/368 * (JDate - 81));
{ Solar longitude }
Slon = 2*PI/24 * (GMT - 12);
{ Solar position vector }
Su = EarthRad*sin(Slon)*cos(Slat);
Sv = EarthRad*cos(Slon)*cos(Slat);
Sw = EarthRad*sin(Slat);
{ Solar cosine }
Scos = (Su*GPu + Sv*GPv + Sw*GPw)/(EarthRad*EarthRad);
{ Earth radiance multiplier }
Eradmult = EarthIrrad / PI * noneg(Scos);
{ cloud cover function }
softif(b, on, off) : if(b-.1, on, if(-.1-b, off, linterp((b+.1)/.2, off, on)));
CloudCover(rho) = if(Cloudy-.02,
softif(Cloudy-(.5+.5*fnoise3(GPu/25,GPv/25,GPw/25)),
CloudRho, rho),
rho);
{ Earth radiance above atmos. }
Eval(rho) = Eradmult * CloudCover(rho);
|