Mengapa perhitungan warna langit saya di Mathematica salah?

17

Saya mencoba menerapkan algoritma untuk menghitung warna langit berdasarkan makalah ini (model Perez '). Sebelum saya mulai memprogram shader saya ingin menguji konsep dalam Mathematica. Sudah ada beberapa masalah yang tidak bisa saya singkirkan. Mungkin seseorang sudah menerapkan algoritma.

Saya mulai dengan persamaan untuk pencahayaan zenital absolut Yz, xzdan yzseperti yang diusulkan dalam makalah (halaman 22). Nilai-nilai untuk Yztampaknya masuk akal. Diagram berikut menunjukkan Yzfungsi jarak zenital matahari untuk kekeruhan T5:

Yz (z)

Fungsi gamma (zenith, azimuth, solarzenith, solarazimuth) menghitung sudut antara titik dengan jarak zenital yang diberikan dan azimuth dan matahari pada posisi yang diberikan. Fungsi ini tampaknya berfungsi juga. Diagram berikut menunjukkan sudut ini untuk solarzenith=0.5dan solarazimuth=0. zenithtumbuh dari atas ke bawah (0 ke Pi / 2), azimuthtumbuh dari kiri ke kanan (-Pi ke Pi). Anda dapat dengan jelas melihat posisi matahari (titik terang, sudut menjadi nol):

gamma (zenith, azimuth, 0,5,0)

Fungsi Perez (F) dan koefisien telah diimplementasikan seperti yang diberikan dalam makalah ini. Maka nilai warna Yxy seharusnya absolute value * F(z, gamma) / F(0, solarzenith). Saya berharap nilai-nilai itu berada dalam kisaran [0,1]. Namun, ini bukan kasus untuk komponen Y (lihat pembaruan di bawah untuk detailnya). Berikut ini beberapa nilai sampel:

{Y, x, y}
{19.1548, 0.25984, 0.270379}
{10.1932, 0.248629, 0.267739]
{20.0393, 0.268119, 0.280024}

Inilah hasil saat ini:

Gambar RGB

Notebook Mathematica dengan semua perhitungan dapat ditemukan di sini dan versi PDF di sini .

Adakah yang tahu apa yang harus saya ubah untuk mendapatkan hasil yang sama seperti di koran?

C suka kode

// 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();
    }
}

Larutan

Seperti yang dijanjikan, saya menulis artikel blog tentang rendering langit. Anda dapat menemukannya di sini .

Nico Schertler
sumber
Saya menduga bahwa lebih banyak orang di sini akan dapat membantu Anda jika Anda mencoba menerapkan algoritma dalam kode aktual (shader atau lainnya) daripada di Mathematica.
Tetrad
2
Ada Mathematica SE , meskipun Anda harus memeriksa FAQ mereka untuk melihat apakah pertanyaan Anda ada pada topik di sana.
MichaelHouse
Nah, pertanyaannya bukan tentang Mathematica, tetapi tentang algoritme. Saya menambahkan versi PDF notebook, sehingga semua orang dapat membacanya. Saya yakin bahwa sintaksinya dapat dipahami oleh programmer umum dan mungkin lebih mudah dipahami daripada kode shader.
Nico Schertler
@NicoSchertler: Masalahnya adalah saya tidak berpikir banyak orang di sini mengerti sintaksis Mathematica. Anda mungkin akan lebih beruntung jika menulis ulang kode Anda dalam bahasa mirip-C atau mirip-Python, setidaknya untuk keperluan pertanyaan ini.
Panda Pajama
2
Pertanyaannya terlalu terlokalisasi dan mungkin ditutup, tetapi terima kasih atas tautannya, ini menarik.
sam hocevar

Jawaban:

4

Ada dua kesalahan dalam matriks yang digunakan untuk xz: 1,00166 harus 0,00166, dan 0,6052 harus 0,06052.

sam hocevar
sumber
Terima kasih atas koreksinya. Hasilnya sekarang terlihat lebih baik, tetapi tidak bisa benar. Saya akan sangat menghargai jika Anda mempertimbangkan pertanyaan yang diperbarui.
Nico Schertler
-2

Sepertinya itu mungkin masalah penskalaan nilai warna?

boobami
sumber
2
Meskipun asumsi Anda mungkin benar, akan lebih bermanfaat untuk memberikan penjelasan tambahan. Karena Anda tidak dapat menjawab seluruh pertanyaan, apa yang Anda tulis harus menjadi komentar di bawah pertanyaan.
danijar
Ini tidak memberikan jawaban untuk pertanyaan itu. Untuk mengkritik atau meminta klarifikasi dari penulis, tinggalkan komentar di bawah posting mereka - Anda selalu dapat mengomentari posting Anda sendiri, dan setelah Anda memiliki reputasi yang cukup Anda akan dapat mengomentari posting apa pun .
MichaelHouse
1
Saya tidak mengerti mengapa saran tidak ditoleransi di sini sama sekali. Jika Anda melihat solusi di atas, itu adalah masalah nilai. Mengarahkan orang ke arah yang benar adalah cara yang lebih baik untuk belajar daripada memberikan solusi tepat setiap saat, bukan? Dan tidak, saya tidak dapat berkomentar di bawah pertanyaannya karena saya tidak diizinkan. Itulah sebabnya saya berkomentar di sini. Tapi terima kasih atas demote-nya. Sangat baik dari Anda dan sangat menyemangati orang baru seperti saya. Terima kasih.
boobami