Почему мой расчет цвета неба в Mathematica неверен?

17

Я пытаюсь реализовать алгоритм для расчета цвета неба на основе этой статьи (модель Переса). Перед тем, как начать программировать шейдер, я хотел проверить концепцию в Mathematica. Уже есть некоторые проблемы, от которых я не могу избавиться. Может быть, кто-то уже реализовал алгоритм.

Я начал с уравнений для абсолютных зенитных яркостей Yz, xzи yzкак было предложено в статье (стр. 22). Значения для Yzкажутся разумными. Следующая диаграмма показывает Yzкак функцию зенитного расстояния солнца для мутности T5:

Yz (г)

Функция gamma (зенит, азимут, соляризенит, соларазимут) вычисляет угол между точкой с заданным зенитным расстоянием и азимутом и солнцем в данной позиции. Эта функция, похоже, тоже работает. Следующая диаграмма показывает этот угол для solarzenith=0.5и solarazimuth=0. zenithрастет сверху вниз (от 0 до Pi / 2), azimuthрастет слева направо (от -Pi до Pi). Вы можете четко видеть положение солнца (яркое пятно, угол становится нулевым):

гамма (зенит, азимут, 0.5,0)

Функция Переса (F) и коэффициенты были реализованы, как указано в статье. Тогда значения цвета Yxy должны быть absolute value * F(z, gamma) / F(0, solarzenith). Я ожидаю, что эти значения будут в диапазоне [0,1]. Однако это не относится к компоненту Y (подробности см. Ниже). Вот несколько примеров значений:

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

Вот текущий результат:

RGB изображение

Блокнот Mathematica со всеми расчетами можно найти здесь, а PDF-версию - здесь .

У кого-нибудь есть идея, что я должен изменить, чтобы получить те же результаты, что и в статье?

C как код

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

Решение

Как и обещал я написал в блоге статью про рендеринг неба. Вы можете найти это здесь .

Нико Шертлер
источник
Я подозреваю, что больше людей здесь смогут помочь вам, если вы попытаетесь реализовать алгоритм в реальном коде (шейдере или иным образом), а не в Mathematica.
Тетрад
2
Есть Mathematica SE , хотя вам нужно будет проверить их ЧаВо, чтобы увидеть, есть ли там ваш вопрос по теме.
MichaelHouse
Ну, вопрос на самом деле не в Mathematica, а в алгоритме. Я добавил PDF-версию записной книжки, чтобы каждый мог ее прочитать. Я уверен, что синтаксис понятен обычному программисту и, вероятно, более понятен, чем шейдерный код.
Нико Шертлер
@NicoSchertler: Проблема в том, что я не думаю, что многие здесь понимают синтаксис Mathematica. Вероятно, вам повезет больше, если вы переписываете свой код на C-подобном или Python-подобном языке, по крайней мере, для целей этого вопроса.
Панда Пижама
2
Вопрос действительно слишком локализован и может быть закрыт, но спасибо за бумажную ссылку, это интересно.
Сэм Хоцевар

Ответы:

4

Есть две ошибки в матрице, используемой для xz: 1.00166 должно быть 0,00166, и 0,6052 должно быть 0,06052.

Сэм Хоцевар
источник
Спасибо за исправление. Результат теперь выглядит лучше, но не может быть правильным. Буду признателен, если вы рассмотрите обновленный вопрос.
Нико Шертлер
-2

Кажется, это может быть проблема масштабирования значения цвета?

boobami
источник
2
Хотя ваше предположение может быть правильным, было бы более полезно дать дополнительные объяснения. Поскольку вы не можете ответить на весь вопрос, то, что вы написали, должно быть комментарием под вопросом.
Данияр
Это не дает ответа на вопрос. Чтобы критиковать или запрашивать разъяснения у автора, оставьте комментарий под его постом - вы всегда можете комментировать свои собственные посты, и, когда у вас будет достаточно репутации, вы сможете комментировать любые посты .
MichaelHouse
1
Я не понимаю, почему предложения здесь вообще не принимаются. Если вы посмотрите на решение выше, это вопрос стоимости. Направлять людей в правильном направлении - лучший способ учиться, чем постоянно давать точные решения, не так ли? И нет, я не могу комментировать ниже его вопрос, потому что мне не разрешено. Вот почему я прокомментировал здесь. Но спасибо за понижение. Очень мило с вашей стороны и очень воодушевляет новых людей, таких как я. Спасибо.
Boobami