Построение векторов, перпендикулярных данному и Ламбертовское распределение
Еще некоторые методы, полезные (если не сказать: необходимые) для рейтрейсеров и других пространственных задачек.
Построение базиса в плоскости с заданной нормалью
Дан единичный вектор n (x,y,z). Нужно получить два единичных вектора np (xp,yp,zp) и nq (xq,yq,zq), таких, что n⊥np,n⊥nq,np⊥nq. Можно также потребовать, чтобы n,np,nq образовывали правый базис.
Вообще говоря, это детерминированная задача, хоть и имеющая множество решений. Но она используется во многих случаях, когда нужно получить вектор, плотность распределения которого обладает осевой симметрией, например, излучить с поверхности «фотон» по Ламбертовскому закону.
Решение.
Сначала находим какой-нибудь вектор t (xt,yt,zt), не коллинеарный с n, например, так:
если |x|>0.58, то t=(0,1,0), иначе t=(1,0,0).
Чаще в литературе встречается условие x≠0, но оно численно неустойчиво - получив на входе вектор n (0.000000001,1,0) , мы можем значительно потерять точность или вообще на каком-то этапе прийти к нулевому вектору и вывалиться с ошибкой при попытке нормализации! Число 0.58 взято почти с потолка, это около 1/√3, и если |x|<0.58, это значит, что одна из оставшихся компонент (или даже обе) заведомо больше, чем |x|. На всякий случай (вдруг эта бОльшая компонента - Y), мы тогда берем вектор (1,0,0). Есть еще параноики, которые находят наименьший коэффициент из трех и ставят единицу именно туда, но, по-моему, это перебор.
Далее все достаточно тривиально:
Метод выглядит неаккуратным, с if'ами и «магическими числами» (0,1,0) и (1,0,0), кажется, что должно существовать решение, которое вообще не «опускается» до координат и целиком выражается в векторной записи, но это не так, и виной всему теорема о причесывании ежа.
Если мы из начала координат построим всевозможные единичные векторы, их концы образуют сферу. К каждому из них мы хотим построить перпендикулярный вектор, и эти векторы, проведенные в каждой точке, образуют касательное векторное поле, причем мы требуем, чтобы оно не обращалось в ноль (в каждой точке было бы направление). Теорема о ежике и гласит, что непрерывным это поле быть не может! И мы вынуждены нарушить изотропность сферы, введя направление просто "наобум", без разницы какое, но главное, чтоб оно было. (1;0;0) и (0;1;0) - достойные претенденты, впрочем, (0;0;1) тоже неплох. Кстати говоря, наиболее быстро построение векторов, перпендикулярных данному будет производиться, если первое векторное умножение выполнить "ручками", с учетом, что две компоненты нулевые.
Излучение по Ламбертовскому закону
У нас есть поверхность с нормалью n (x,y,z). Идеально матовая поверхность отражает свет во все стороны независимо от того, откуда этот свет пришел, причем интенсивность излучения в направлении γ от нормали пропорциональна cosγ, такое распределение по полусфере и называется Ламбертовским (Lambertian). В рейтрейсерах, работающих по методу Монте-Карло, вместо излучения во все стороны случайным образом выбирается одно направление, куда излучить «фотон», при этом плотность распределения по сфере также должна быть Ламбертовской.
Эта задача проще решается в сферических координатах, только в качестве оси Y выступает вектор n, а осей X,Z - np,nq, полученные по алгоритму, приведенному выше. Излучение обладает осевой симметрией, поэтому азимут θ является равномерно распределенной на [0;2π) случайной величиной, т.е θ=2πη.
Найдем плотность распределения для угла места φ. Рассмотрим на сфере полосу, ограниченную окружностями φ, φ+dφ, вероятность, что вектор направления попадет в эту полосу, пропорциональна площади полосы, помноженной на sinφ, т.е
p[φ≤ξ<φ+dφ]=k*cosφdφ*sinφ,
где k-коэффициент пропорциональности. Плотность распределения равна, очевидно:
fφ(x)=k*cosxsinx
Найдем функцию распределения:
Возникающие по пути константы мы включаем в коэффициент пропорциональности, от того он все время разный. Отнормируем его из условия Fφ(π/2)=1, получим k1=1 и Fφ(x)=sin2x
Если у нас есть равномерно распределенная на [0;1) случайная величина ξ, то угол места с распределением Fφ(x) выражается через нее так:
В декартовых координатах (в базисе из векторов n,np,nq) получим:
Пару (cos(2πη),sin(2πη)) можно сформировать по алгоритму, описанному в разделе «выбор точки на окружности», что позволит обойтись без тригонометрических функций. С другой стороны, в FPU x86 есть команда FSINCOS, которая считает синус и косинус от числа одним махом, так что однозначно указать самый быстрый метод не представляется возможным.
Наконец, запишем вектор направления излучения:
nизл=xnp+yn+znq
Похожим образом можно генерировать другие осесимметричные распределения, например, для моделирования глянцевой поверхности можно получить сначала вектор зеркального отражения, а потом, построив вокруг него базис, использовать распределения Фонга cosnγ, с n>1 (n=1 возвращает нас к Ламберту).
Построение базиса в плоскости с заданной нормалью
Дан единичный вектор n (x,y,z). Нужно получить два единичных вектора np (xp,yp,zp) и nq (xq,yq,zq), таких, что n⊥np,n⊥nq,np⊥nq. Можно также потребовать, чтобы n,np,nq образовывали правый базис.
Вообще говоря, это детерминированная задача, хоть и имеющая множество решений. Но она используется во многих случаях, когда нужно получить вектор, плотность распределения которого обладает осевой симметрией, например, излучить с поверхности «фотон» по Ламбертовскому закону.
Решение.
Сначала находим какой-нибудь вектор t (xt,yt,zt), не коллинеарный с n, например, так:
если |x|>0.58, то t=(0,1,0), иначе t=(1,0,0).
Чаще в литературе встречается условие x≠0, но оно численно неустойчиво - получив на входе вектор n (0.000000001,1,0) , мы можем значительно потерять точность или вообще на каком-то этапе прийти к нулевому вектору и вывалиться с ошибкой при попытке нормализации! Число 0.58 взято почти с потолка, это около 1/√3, и если |x|<0.58, это значит, что одна из оставшихся компонент (или даже обе) заведомо больше, чем |x|. На всякий случай (вдруг эта бОльшая компонента - Y), мы тогда берем вектор (1,0,0). Есть еще параноики, которые находят наименьший коэффициент из трех и ставят единицу именно туда, но, по-моему, это перебор.
Далее все достаточно тривиально:
Метод выглядит неаккуратным, с if'ами и «магическими числами» (0,1,0) и (1,0,0), кажется, что должно существовать решение, которое вообще не «опускается» до координат и целиком выражается в векторной записи, но это не так, и виной всему теорема о причесывании ежа.
Если мы из начала координат построим всевозможные единичные векторы, их концы образуют сферу. К каждому из них мы хотим построить перпендикулярный вектор, и эти векторы, проведенные в каждой точке, образуют касательное векторное поле, причем мы требуем, чтобы оно не обращалось в ноль (в каждой точке было бы направление). Теорема о ежике и гласит, что непрерывным это поле быть не может! И мы вынуждены нарушить изотропность сферы, введя направление просто "наобум", без разницы какое, но главное, чтоб оно было. (1;0;0) и (0;1;0) - достойные претенденты, впрочем, (0;0;1) тоже неплох. Кстати говоря, наиболее быстро построение векторов, перпендикулярных данному будет производиться, если первое векторное умножение выполнить "ручками", с учетом, что две компоненты нулевые.
Излучение по Ламбертовскому закону
У нас есть поверхность с нормалью n (x,y,z). Идеально матовая поверхность отражает свет во все стороны независимо от того, откуда этот свет пришел, причем интенсивность излучения в направлении γ от нормали пропорциональна cosγ, такое распределение по полусфере и называется Ламбертовским (Lambertian). В рейтрейсерах, работающих по методу Монте-Карло, вместо излучения во все стороны случайным образом выбирается одно направление, куда излучить «фотон», при этом плотность распределения по сфере также должна быть Ламбертовской.
Эта задача проще решается в сферических координатах, только в качестве оси Y выступает вектор n, а осей X,Z - np,nq, полученные по алгоритму, приведенному выше. Излучение обладает осевой симметрией, поэтому азимут θ является равномерно распределенной на [0;2π) случайной величиной, т.е θ=2πη.
Найдем плотность распределения для угла места φ. Рассмотрим на сфере полосу, ограниченную окружностями φ, φ+dφ, вероятность, что вектор направления попадет в эту полосу, пропорциональна площади полосы, помноженной на sinφ, т.е
p[φ≤ξ<φ+dφ]=k*cosφdφ*sinφ,
где k-коэффициент пропорциональности. Плотность распределения равна, очевидно:
fφ(x)=k*cosxsinx
Найдем функцию распределения:
Возникающие по пути константы мы включаем в коэффициент пропорциональности, от того он все время разный. Отнормируем его из условия Fφ(π/2)=1, получим k1=1 и Fφ(x)=sin2x
Если у нас есть равномерно распределенная на [0;1) случайная величина ξ, то угол места с распределением Fφ(x) выражается через нее так:
В декартовых координатах (в базисе из векторов n,np,nq) получим:
Пару (cos(2πη),sin(2πη)) можно сформировать по алгоритму, описанному в разделе «выбор точки на окружности», что позволит обойтись без тригонометрических функций. С другой стороны, в FPU x86 есть команда FSINCOS, которая считает синус и косинус от числа одним махом, так что однозначно указать самый быстрый метод не представляется возможным.
Наконец, запишем вектор направления излучения:
nизл=xnp+yn+znq
Похожим образом можно генерировать другие осесимметричные распределения, например, для моделирования глянцевой поверхности можно получить сначала вектор зеркального отражения, а потом, построив вокруг него базис, использовать распределения Фонга cosnγ, с n>1 (n=1 возвращает нас к Ламберту).