Interpolate Gradient noise в Julia

Question

Я пытаюсь реализовать и изучить код с here, но у меня проблемы с интерполяцией. Это мой код из моей реализации в Джулии:

lerp(a, b, w) = a * (1 - w) + b * w 

function Noise(x, y) 
    n = x + y * 53; 
    n = (n << 13) $ n; 
    return (1.0 - ((n * ((n * n * 15731) + 789221) + 1376312589) & 0x7fffffff)/1073741824.0); 
end 

function coherentNoise(x,y) 
    x₁ = convert(Int64,modf(x)[2]) 
    y₁ = convert(Int64,modf(y)[2]) 
    xᵣ = x - x₁ 
    yᵣ = y - y₁ 
    q1 = Noise(x₁-1,y₁+1) 
    q2 = Noise(x₁+1,y₁+1) 
    q3 = Noise(x₁-1,y₁-1) 
    q4 = Noise(x₁-1,y₁+1) 

    v = lerp(q1,q2,xᵣ) 
    v1 = lerp(q3,q4,xᵣ) 
    return abs(lerp(v,v1,yᵣ)) 
    #value = bilinear(q1,q2,q3,q4,x₁,x₁-1,x₁+1,y₁,y₁-1,y₁+1) 
    #return abs(value) 
end 


arr = Array{Float64}(height,width) 

for x = 1:height 
    for y = 1:width 
     nx = x/60 
     ny = y/60 
     arr[x,y]=coherentNoise(nx,ny) 
    end 
end 

это результат:

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

# bilinear interpolation 
function bilinear(Q1,Q2,Q3,Q4,x,x1,x2,y,y1,y2) 
    R1 = ((x2-x)/(x2-x1))*Q3 + ((x-x1)/(x2-x1))*Q4 
    R2 = ((x2-x)/(x2-x1))*Q1 + ((x-x1)/(x2-x1))*Q2 
    P = ((y2-y)/(y2-y1))*R1 + ((y-y1)/(y2-y1))*R2 
    return P 
end 

результат:

Я не знаю, если это проблема с Noise функции или что-то другое.

Update

Я получения наилучших результатов с помощью билинейной интерполяции со следующей функцией генератора шума:

function coherentNoise(x,y,n::Int64) 
    #octa3e 
    o = 0.25 
    #octave x and y 
    x /=n 
    y /=n 
    #decimal value 
    x₁ = x < 0 ? x : x-o 
    y₁ = y < 0 ? y : y-o 
    x₂ = x+o 
    y₂ = y+o 
    #calculate corners data 
    q1 = abs(Noise(floor(Int64,x₁),floor(Int64,y₂))) 
    q2 = abs(Noise(floor(Int64,x₂),floor(Int64,y₂))) 
    q3 = abs(Noise(floor(Int64,x₁),floor(Int64,y₁))) 
    q4 = abs(Noise(floor(Int64,x₂),floor(Int64,y₁))) 
    #bilenear interpolation 

    value = bilinear(q2,q1,q4,q3,x,x₁,x₂,y,y₁,y₂) 
    return abs(value) 
end 

Результат:

Я думаю, что проблема как я выбираю точки для билинейной интерполяции, не знаю, как это сделать правильно.

Update 2:

Это код для создания изображений:

for x = 1:height 
    for y = 1:width 
     arr[x,y]=coherentNoise(x,y,50) 
    end 
end 

imwrite(convert(Image,arr),"desktop/projects/julia/Noise/test2.png") 

Update 3

Благодаря ответу Dan Getz я получаю лучший выход, но с некоторыми странно "черви"

Answer 1

Без возможности проверки кода это должно быть подтверждено. В принципе, я думаю, что ваша интуиция является проблемой в выборе точек интерполяции. Это попытка исправить это (см. Также комментарии).

function coherentNoise(x,y) 
    x₁,xᵣ = floor(Int64,x),mod(x,1) 
    y₁,yᵣ = floor(Int64,y),mod(y,1) 
    q00 = Noise(x₁,y₁) 
    q10 = Noise(x₁+1,y₁) 
    q01 = Noise(x₁,y₁+1) 
    q11 = Noise(x₁+1,y₁+1) 

    v0 = lerp(q00,q10,xᵣ) 
    v1 = lerp(q01,q11,xᵣ) 
    noise = lerp(v0,v1,yᵣ) 
    return noise 
end 

Примечания:

• q ы переименован в q00 и тому подобное, чтобы указать, какой угол от смещения Х (первое число) и Y (второй номер).
• Исправлены идентичные q1 и q4.
• floor и mod(,1) используются вместо modf, так как они выбирают согласованное направление округления (вниз) независимо от знака.
• Функция Noise должна указывать аннотацию ее параметров, чтобы ограничить значения Int32, а ; на линии заканчивается этой функцией. Это хеш-функция, и тот же эффект может быть легко закодирован с помощью функции хеши Джулии и небольшой оберткой.

Например:

Noise(x,y) = 1.0-2.0*hash((x,y))/typemax(UInt) 

является более Джулией стиля капли в замене Noise.

Если вы попробуете это и прокомментируете ошибки, мы можем попытаться выполнить функцию fixup.