Ich versuche einen Algorithmus zu implementieren, um die Himmelsfarbe basierend auf diesem Papier (Perez 'Modell) zu berechnen . Bevor ich den Shader programmiere, wollte ich das Konzept in Mathematica testen. Es gibt bereits einige Probleme, die ich nicht loswerden kann. Vielleicht hat jemand den Algorithmus schon implementiert.
Ich begann mit Gleichungen für den absoluten zenital Helligkeiten Yz
, xz
und yz
wie in dem Papier (Seite 22) vorgeschlagen. Die Werte für Yz
scheinen vernünftig zu sein. Das folgende Diagramm zeigt Yz
in Abhängigkeit vom Zenitabstand der Sonne bei einer Trübung T
von 5:
Die Funktion Gamma (Zenit, Azimut, Solarzenit, Solarazimut) berechnet den Winkel zwischen einem Punkt mit dem angegebenen Zenitabstand und Azimut und der Sonne an der angegebenen Position. Diese Funktion scheint auch zu funktionieren. Das folgende Diagramm zeigt diesen Winkel für solarzenith=0.5
und solarazimuth=0
. zenith
wächst von oben nach unten (0 bis Pi / 2), azimuth
wächst von links nach rechts (-Pi bis Pi). Sie können den Sonnenstand deutlich sehen (der helle Punkt, der Winkel wird zu Null):
Die Perez-Funktion (F) und die Koeffizienten wurden wie im Artikel angegeben implementiert. Dann sollten die Farbwerte Yxy sein absolute value * F(z, gamma) / F(0, solarzenith)
. Ich erwarte, dass diese Werte im Bereich [0,1] liegen. Dies ist jedoch bei der Y-Komponente nicht der Fall (Details siehe Update unten). Hier sind einige Beispielwerte:
{Y, x, y}
{19.1548, 0.25984, 0.270379}
{10.1932, 0.248629, 0.267739]
{20.0393, 0.268119, 0.280024}
Hier ist das aktuelle Ergebnis:
Das Mathematica-Notizbuch mit allen Berechnungen finden Sie hier und die PDF-Version hier .
Hat jemand eine Idee, was ich ändern muss, um die gleichen Ergebnisse wie in der Zeitung zu erzielen?
C wie Code
// this function returns the zenital Y component for
// a given solar zenital distance z and turbidity T
float Yz(float z, float T)
{
return (4.0453 * T - 4.9710)*tan( (4.0f/9-T/120)*(Pi-2*z) ) - 0.2155 * T + 2.4192
}
// returns zenital x component
float xz(float z, float T)
{
return //matrix calculation, see paper
}
// returns zenital y component
float yz(float z, float T)
{
return //matrix calculation, see paper
}
// returns the rgb color of a Yxy color
Color RGB(float Y, float x, float y)
{
Matrix m; //this is a CIE XYZ -> RGB conversion matrix
Vector v;
v.x = x/y*Y;
v.y = Y;
v.z = (1-x-y)/y*Y;
v = M * v; //matrix-vector multiplication;
return Color ( v.x, v.y, v.z );
}
// returns the 5 coefficients (A-E) for the given turbidity T
float[5] CoeffY(float T)
{
float[5] result;
result[0] = 0.1787 * T - 1.4630;
result[1] = -0.3554 * T + 0.4275;
...
return result;
}
//same for Coeffx and Coeffy
// returns the angle between an observed point and the sun
float PerezGamma(float zenith, float azimuth, float solarzenith, float solarazimuth)
{
return acos(sin(solarzenith)*sin(zenith)*cos(azimuth-solarazimuth)+cos(solarzenith)*cos(zenith));
}
// evalutes Perez' function F
// the last parameter is a function
float Perez(float zenith, float gamma, float T, t->float[5] coeffs)
{
return (1+coeffs(T)[0] * exp(coeffs(T)[1]/cos(zenith)) *
(1+coeffs(T)[2] * exp(coeffs(T)[3]*gamma) +
coeffs(T)[4]*pow(cos(gamma),2))
}
// calculates the color for a given point
YxyColor calculateColor(float zenith, float azimuth, float solarzenith, float solarazimuth, float T)
{
YxyColor c;
float gamma = PerezGamma(zenith, azimuth, solarzenith, solarazimuth);
c.Y = Yz(solarzenith, T) * Perez(zenith, gamma, T, CoeffY) / Perez(0, solarzenith, T, CoeffY);
c.x = xz(solarzenith, T) * Perez(zenith, gamma, T, Coeffx) / Perez(0, solarzenith, T, Coeffx);
c.y = yz(solarzenith, T) * Perez(zenith, gamma, T, Coeffy) / Perez(0, solarzenith, T, Coeffy);
return c;
}
// draws an image of the sky
void DrawImage()
{
for(float z from 0 to Pi/2) //zenithal distance
{
for(float a from -Pi to Pi) //azimuth
{
YxyColor c = calculateColor(zenith, azimuth, 1, 0, 5);
Color rgb = RGB(c.Y, c.x, c.y);
setNextColor(rgb);
}
newline();
}
}
Lösung
Wie versprochen habe ich einen Blogartikel über das Rendern des Himmels geschrieben. Sie finden es hier .
Antworten:
Die verwendete Matrix enthält zwei Fehler
xz
: 1.00166 sollte 0.00166 sein und 0.6052 sollte 0.06052 sein.quelle
Es scheint, als könnte es sich um ein Problem mit der Farbwertskalierung handeln.
quelle