Бесплатная реклама

Автор: Ильин Михаил «Mikle»

Опубликовано: 01.01.1970

Изменено: 01.01.1970

Постоянная ссылка

Комментарии [2]

Имитация раскачивания дерева под воздействием ветра


Статья демонстрирует простой метод имитации воздействия ветра на дерево, реализованный на вертексном шейдере.


Содержание


Шаг 1. Подготовка скелета дерева
Шаг 2. Имитация ветра
Шаг 3. Добавляем хаотичность и инерцию
Демо-версия и исходный код


Статья демонстрирует простой метод имитации воздействия ветра на дерево (куст), реализованный на вертексном шейдере. К статье прилагаются три примера, содержащие скомпилированные EXE файлы и их исходники на языке Visual Basic из состава Microsoft Visual Studio 2008 Express. Для запуска EXE файлов необходимо иметь установленными NET Framework 2.0, DirectX9.0c и Managed DirectX. Видеоадаптер должен поддерживать вертексные шейдеры vs_1_1.

Демонстрационные примеры следуют один за другим, шаг за шагом описывая построение программы.


Шаг 1. Подготовка скелета дерева


Это подготовительный этап, где мы создаем скелет дерева, с которым будем работать в дальнейшем. То есть мы будем использовать не готовую модель, а только ее имитацию, достаточную для оценки работы физики. Дерево можно создать на основе фрактала по такой схеме:
  • Имеем 2 вектора, начальную позицию и направление.
  • Позиция задает начало, а сумма позиции и направления конец отрезка, который будет добавлен в нашу модель.
  • Считаем новой позицией конец отрезка.
  • Трижды вычисляем новое направление, отклоняя старое направление в разные стороны на некоторый угол, длина вектора нового направления должна быть несколько короче старого (вводим коэффициент). Рекурсивно переходим к п.2 для каждой новой пары позиция-направление.

Таким образом, первый отрезок становится стволом дерева, от ствола отходит три ветви, несколько короче самого ствола, от каждой ветви еще три ветви и т. д. Глубину рекурсии установим в 7.

На языке Visual Basic это выглядит так:

Private Sub AddStick(ByVal Level As Integer, ByVal pV As Integer, ByVal Dir As Vector3)
    Dim d, vS, vC As Vector3
    Dim p As Integer, s As Single
    pVert = pVert + 1
    Vert(pVert).Position = Vert(pV).Position + Dir
    Ind(pInd) = pV
    pInd = pInd + 1
    Ind(pInd) = pVert
    pInd = pInd + 1

    If Level = 1 Then Exit Sub

    StickLen = StickLen * 0.7
    d = Vector3.Normalize(Dir) * StickLen
    s = StickLen * 0.6
    vS = Vector3.Normalize(Vector3.Cross(Dir, New Vector3(Rnd, Rnd, Rnd))) * s
    vC = Vector3.Normalize(Vector3.Cross(Dir, vS)) * s

    p = pVert
    AddStick(Level - 1, p, vS * 0.5 + vC * 0.866 + d * (Rnd() + 1.5) * 0.5)
    AddStick(Level - 1, p, vS * -1 + d * (Rnd() + 1.5) * 0.5)
    AddStick(Level - 1, p, vS * 0.5 + vC * -0.866 + d * (Rnd() + 1.5) * 0.5)
    StickLen = StickLen / 0.7
End Sub


Аргументами процедуры являются максимальный уровень рекурсии, указатель на вертекс, содержащий координаты позиции и направление. Для того, чтобы дерево не было фракталом идеальной формы, а выглядело естественно, в длину и направление вычисляемых векторов вносятся небольшие случайные поправки.

Более подробно останавливаться на генерации фрактала не буду, тема статьи другая, а в Интернете достаточно информации по этому вопросу.

Формат вертекса модели пока содержит только Position, для отрисовки используется простейший вертексный шейдер:

  vs_1_1

  dcl_position v0
  def c8, 1, 1, 1, 0

  m4x4 oPos, v0, c0
  mov oD0, c8


Просто перерасчитываем координаты и выводим в oD0 белый цвет, заданный в константе c8. В константы c0-c3 загружаем матрицу, обобщающую трансформации проекции, видовую и мировую:

    Dim M1 As Matrix
    Dim M2 As Matrix

    M1 = Dev.Transform.World
    M2 = Dev.Transform.View
    M1 = Matrix.Multiply(M1, M2)
    M2 = Dev.Transform.Projection
    M1 = Matrix.Multiply(M1, M2)
    M2 = Matrix.TransposeMatrix(M1)
    Dev.SetVertexShaderConstant(0, M2)


Вот и всё, посмотрите пример в папке Step1, в нем двойным щелчком мыши можно генерировать новое дерево, левой кнопкой вращаем камеру вокруг дерева, которое пока статично. Уровень рекурсии 7, что соответствует 1094 вертексов (1093 ветви).

http://www.uraldev.ru/files/articles/tree-shaking-simulation-by-wind/fractal_tree.jpg




Шаг 2. Имитация ветра


Теперь немного теории. Зададимся целью вносить в координаты вертексов некоторую регулируемую поправку, имитирующую сгибание ветки от ветра. Можно вычислять матрицу поворота вокруг оси, проходящей через начало ветки и перпендикулярной вектору силы (ветра), это было бы правильно, но очень громоздко. При уровне рекурсии 7 крайние ветви должны испытывать отклонение по семи матрицам! Это тяжелые расчеты, и это сильно «потяжелевший» вертекс.

Попробуем найти более простой, пусть и не столь точный метод. Будем просто смещать вершину ветки в нужном направлении, наша задача найти это направление. Ветер не может быть вертикальным, то есть наша задача найти вектор смещения для ветра единичной силы для каждой ветки для двух направлений ветра, - вдоль оси X и вдоль оси Z, - и записать в вертексы эти вектора. Далее в шейдере останется умножить первый вектор на компоненту X ветра, второй на компоненты Z и просуммировать полученные вектора - получим результирующее смещение.

Смещение может быть только перпендикулярно направлению ветки (ведь ветка не может изменять длину), оно должно лежать в одной плоскости с векторами направлений ветки и ветра. Чтобы найти такой вектор, достаточно найти Cross продукт от вектора, перпендикулярного ветке и ветру, и вектора самой ветки. А перпендикулярный вектор так же находится Cross продуктом векторов ветра и ветки:

  v = Vector3.Normalize(Vector3.Cross(Dir, Vector3.Cross(New Vector3(1, 0, 0), Dir))) * Dir.Length


Длина вектора смещения должна быть пропорциональна длине ветки, поэтому мы нормализовали результат и умножили на длину ветки.

Это для ветра вдоль оси X. То же самое для ветра вдоль оси Z:

  v = Vector3.Normalize(Vector3.Cross(Dir, Vector3.Cross(New Vector3(0, 0, 1), Dir))) * Dir.Length


Чтобы учитывать то, что кроме собственных колебаний, ветка подвержена колебаниям несущих ветвей, будем перед записью в вертекс вектора суммировать его с соответствующим вектором несущей ветки, а тот вектор уже просуммирован с соответствующим вектором своей несущей ветки и т. д. до основания ствола, которое, естественно, неподвижно. Кроме того, введем зависимость от толщины ветки – чем ветка тоньше, тем она легче гнется, будем умножать вектора не на длину, а на некоторую величину, зависящую от уровня ветки в рекурсии, можно задать такую величину таблично, а можно и подобрать формулу имперически, например так:

  s = Dir.Length * (0.3 + meLevels - Level)


Здесь meLevels –  уровень рекурсии ствола, а Levels – уровень рекурсии текущей ветки.

Теперь в вертексе кроме координат содержатся еще два вектора, задающие смещения при ветре вдоль осей X и Z. В вертексный шейдер их можно передать, например, как нормаль и бинормаль (а можно как-нибудь по-другому, вектор – он и есть вектор). Направление ветра передаем в шейдер через константу c4:

  Dim w As Vector4
  w.X = Wind.X
  w.Y = 0
  w.Z = Wind.Z
  w.W = 0
  Dev.SetVertexShaderConstant(4, w)


А сам шейдер выглядит так:

  vs_1_1

  dcl_position v0
  dcl_normal v1
  dcl_binormal v2

  def c5, 1, 1, 1, 0

  mov r0.w, c5.w
  mul r0.xyz, v1.xyz, c4.x
  add r1, r0, v0

  mul r0.xyz, v2.xyz, c4.z
  add r1, r1, r0

  m4x4 oPos, r1, c0
  mov oD0, c5


Тут вычисляется, как уже было сказано чуть выше, сумма произведений векторов-коэффициентов на соответствующие  компоненты вектора ветра, компонента w приравнивается нулю, и само смещение суммируется с координатами вершины.

В примере Step2, как и в предыдущем, левой кнопкой вращаем камеру, правой направление ветра, а клавишами  1 – 9 задаем силу ветра. Мы видим сгибающееся в разных направлениях дерево, но на ветер это пока непохоже, не хватает хаотичности. При смене силы ветра дерево меняет вид мгновенно – то есть не хватает еще и инерции.


Шаг 3. Добавляем хаотичность и инерцию


Хаотичность (неравномерность) ветра хорошо описывает так называемый шум Перлина. Создадим процедуры для его генерации и интерполяции для выборки промежуточных значений:

  Private Function GetNoise(ByVal s As Single) As Single
    Dim n1, n2 As Integer
    Dim t As Single
    n1 = Int(s)
    t = s - n1
    n1 = n1 And 255
    n2 = (n1 + 1) And 255
    GetNoise = Noise(n1) * (1 - t) + Noise(n2) * t
  End Function

  Private Sub GenNoise()
    Dim n As Integer, dn As Integer
    dn = 32
    Do
      For n = dn To 255 Step dn * 2
        Noise(n) = (Noise(n - dn) + Noise((n + dn) And 255)) * 0.5 + (Rnd() - 0.5) * dn / 16
      Next n
      dn = dn >> 1
      If dn = 0 Then Exit Do
    Loop
  End Sub


Подробнее на этом останавливаться не буду, так как про шумы Перлина в интернете информации достаточно.

Инерция подразумевает постепенное изменение позиции при мгновенно изменившейся силе, для ее непосредственной имитации необходимо результат работы шейдера сохранять в вертексе до следующего кадра, во-первых, это заняло бы много памяти, во-вторых у нас нет такой возможности, поэтому опять будет фейк. Постепенно изменять будем сам ветер. Можно использовать обычную линейную интерполяцию между предыдущим и текущим значениями ветра по времени:

  sv4 = sv4 * (1 - dt) + dv4 * dt


Здесь sv4 – статичный (сохраняющий значение между вызовами процедуры) вектор, dv4 – вектор ветра, а dt – время, прошедшее с предыдущего кадра. И сразу внесем хаотичность, добавим к sv4 случайную величину, заданную шумом Перлина:

  dv4.X = v.X * (1 - 0.5 * GetNoise(t * 63))
  dv4.Y = 0
  dv4.Z = v.Z * (1 - 0.5 * GetNoise(t * 50))
  dv4.W = 0


Константы 63 и 50 подобраны эмпирически, это скорость движения по шуму компонент X и Z, они заданы разными, таким образом, ветер будет менять не только силу, но и направление.

Уже стало реальнее, но добавим последний штрих – резонансные случайные движения тонких веток, зависящие от силы ветра. Добавим в вертекс три коэффициента, вместе они прекрасно уместятся в незанятый еще Tangent. Для трех верхних ступеней рекурсии заполним их случайными числами так, чтобы для тонких веток значения были больше:

    If Level = 1 Then
      s = Dir.Length * 0.3
    ElseIf Level = 2 Then
      s = Dir.Length * 0.15
    ElseIf Level = 3 Then
      s = Dir.Length * 0.075
    Else
      s = 0
    End If
    v = New Vector3(s * (Rnd() - 0.5), s * (Rnd() - 0.5), s * (Rnd() - 0.5))


Мы внесли сразу три коэффициента, будем отправлять в шейдер три соответствующие выборки шума:

    v4.X = GetNoise(QTime() * 33) * 10 * s
    v4.Y = GetNoise(QTime() * 52) * 10 * s
    v4.Z = GetNoise(QTime() * 94) * 10 * s
    Dev.SetVertexShaderConstant(5, v4)


Константы 33, 52, 94 определяют скорость дрожания веток, можно подобрать другие величины на свой вкус. Константа 10 определяет амплитуду колебаний, а параметр s – это сила ветра (длина вектора).

Мы не зря записали в вертексы случайные коэффициенты, на колебания разных веток выборки шума будут сказываться в разной степени, ветки будут двигаться не синхронно.Наша задача в шейдере перемножить соответствующие коэффициенты на выборки шума и просуммировать эти величины.

С этой задачей прекрасно справится единственная инструкция dp3! Итак, наш шейдер:

  vs_1_1

  dcl_position v0
  dcl_normal v1
  dcl_binormal v2
  dcl_tangent v3

  def c6, 1, 1, 1, 0

  dp3 r2, v3, c5

  add r3, r2, c4.xxxw
  mul r0, v1, r3
  add r1, r0, v0

  add r3, r2, c4.zzzw
  mul r0, v2, r3
  add r1, r1, r0

  m4x4 oPos, r1, c0
  mov oD0, c6


В примере Step3 глубина рекурсии увеличена до 8 (3281 веток), управление такое же, как и в Step2. Реализованы инерция, неравномерность и случайные колебания тонких веток.

Статья окончена, в качестве творческого задания предлагаю разработать алгоритм записи таких же коэффициентов в уже готовую, созданную моделлерами, модель дерева.


Демо-версия и исходный код


Скачать демо-версию и исходный код можно здесь: wind_tree_simulation.zip






Статья участвовала в конкурсе статей по программированию #3 (2008).