melhoria da precisão da interseção de raios e elipsóides

Preciso aprimorar a precisão para a função em um dos meusShader de fragmento de dispersão atmosférica GLSL que calcula a interseção entre o raio único e o elipsóide alinhado ao eixo.

Essa é a principal função do shader de dispersão atmosférica da mina. O antigo shader original estava ligadofloats e para renderização normal foi bom, mas após a adição do zoom, descobri que, com distâncias relativamente pequenas, a precisão se perde. Em carros alegóricos, as distâncias utilizáveis para a Terra eram de apenas 0,005 UA (Unidade Astronômica). Então, tentei portar a função crítica paradouble e ajuda, agora a distância utilizável é de cerca de 1,0 AU (com pequenos artefatos)

Isto é odouble versão da função dentro do Fragment Shader (a fonte antiga usa coisas obsoletas !!!)

#extension GL_ARB_gpu_shader_fp64 : enable
double abs(double x) { if (x<0.0) x=-x; return x; }
// compute length of ray(p0,dp) to intersection with ellipsoid((0,0,0),r) -> view_depth_l0,1
// where r.x is elipsoid rx^-2, r.y = ry^-2 and r.z=rz^-2
float view_depth_l0=-1.0,view_depth_l1=-1.0;
bool _view_depth(vec3 _p0,vec3 _dp,vec3 _r)
    {
    double a,b,c,d,l0,l1;
    dvec3 p0,dp,r;
    p0=dvec3(_p0);
    dp=dvec3(_dp);
    r =dvec3(_r );
    view_depth_l0=-1.0;
    view_depth_l1=-1.0;
    a=(dp.x*dp.x*r.x)
     +(dp.y*dp.y*r.y)
     +(dp.z*dp.z*r.z); a*=2.0;
    b=(p0.x*dp.x*r.x)
     +(p0.y*dp.y*r.y)
     +(p0.z*dp.z*r.z); b*=2.0;
    c=(p0.x*p0.x*r.x)
     +(p0.y*p0.y*r.y)
     +(p0.z*p0.z*r.z)-1.0;
    d=((b*b)-(2.0*a*c));
    if (d<0.0) return false;
    d=sqrt(d);
    l0=(-b+d)/a;
    l1=(-b-d)/a;
    if (abs(l0)>abs(l1)) { a=l0; l0=l1; l1=a; }
    if (l0<0.0)          { a=l0; l0=l1; l1=a; }
    if (l0<0.0) return false;
    view_depth_l0=float(l0);
    view_depth_l1=float(l1);
    return true;
    }
entrada é raio e raios ^ -2 de elipsóide

saída é a distância entre p0 e cruzamentos

a precisão das variáveis de entrada e saída é flutuante (é suficiente)

É assim que ele cuida da portabilidade para o Double

Então a pergunta é: Q1. Como melhorar a precisão dessa função?

precisão alvo paraview_depth_l0,view_depth_l1 é+/- 20 m para distância|p0|=100 AU

Isso seria ideal, agora parece que +/- 5 km para uma distância de 10 UA, o que é ruim. Até 10 vezes a computação precisa será um grande passo adiante, alguma idéia?

[edit1] l0, l1 range

Eu estava errado flutuar conversão deview_depth_l0,view_depth_l1 é a causa dos artefatos. Depois de mudar para uma distância relativa, a precisão melhorou bastante. Acabei de adicionar isso:

    // relative shift to preserve accuracy
    const double m0=1000000000.0; // >= max view depth !!!
    if (l0>m0){ a=floor(l0/m0)*m0; a-=m0; if (l1>l0) l1-=a; l0-=a; }

antes disso:

    view_depth_l0=float(l0);
    view_depth_l1=float(l1);
    return true;
    }

O resto da alça do shaderl0,l1 como valores relativos de qualquer maneira, agora o resultado é o seguinte:

para distâncias de até 10,0 UA, está bom agora (artefatos são perceptíveis apenas em zooms muito altos), novos artefatos são causados provavelmente em outros lugares; portanto, será necessário pesquisar mais quando tiver tempo e vontade.

[edit2] deslocando p0 ao longo de dp para mais perto (0,0,0)

A implementação precisa de funções de normalização e comprimento relativamente caras e o resultado sem deslocamento de faixa (edição1) foi um pouco melhor que a função bruta, mas a melhoria não é muito grande. Com a mudança de faixa (edit1), o resultado é o mesmo de antes, portanto não é esse o caminho. Minha conclusão é que todos os artefatos restantes não são causados pela função do departamento de exibição em si.

Vou tentar portar shader para#version 400 + fp64 em tudo para verificar se os dados de entrada não são arredondados por flutuar muito

[edit3] código fonte real

#extension GL_ARB_gpu_shader_fp64 : enable
double abs(double x) { if (x<0.0) x=-x; return x; }
// compute length of ray(p0,dp) to intersection with ellipsoid((0,0,0),r) -> view_depth_l0,1
// where r.x is elipsoid rx^-2, r.y = ry^-2 and r.z=rz^-2
float view_depth_ll= 0.0, // shift to boost accuracy
      view_depth_l0=-1.0, // view_depth_ll+view_depth_l0 first hit
      view_depth_l1=-1.0; // view_depth_ll+view_depth_l1 second hit
const double view_depth_max=100000000.0; // > max view depth
bool _view_depth(vec3 _p0,vec3 _dp,vec3 _r)
    {
    dvec3 p0,dp,r;
    double a,b,c,d,l0,l1;
    view_depth_ll= 0.0;
    view_depth_l0=-1.0;
    view_depth_l1=-1.0;
    // conversion to double
    p0=dvec3(_p0);
    dp=dvec3(_dp);
    r =dvec3(_r );
    // quadratic equation a.l.l+b.l+c=0; l0,l1=?;
    a=(dp.x*dp.x*r.x)
     +(dp.y*dp.y*r.y)
     +(dp.z*dp.z*r.z);
    b=(p0.x*dp.x*r.x)
     +(p0.y*dp.y*r.y)
     +(p0.z*dp.z*r.z); b*=2.0;
    c=(p0.x*p0.x*r.x)
     +(p0.y*p0.y*r.y)
     +(p0.z*p0.z*r.z)-1.0;
    // discriminant d=sqrt(b.b-4.a.c)
    d=((b*b)-(4.0*a*c));
    if (d<0.0) return false;
    d=sqrt(d);
    // standard solution l0,l1=(-b +/- d)/2.a
    a*=2.0;
    l0=(-b+d)/a;
    l1=(-b-d)/a;
    // alternative solution q=-0.5*(b+sign(b).d) l0=q/a; l1=c/q; (should be more accurate sometimes)
//  if (b<0.0) d=-d; d=-0.5*(b+d);
//  l0=d/a;
//  l1=c/d;
    // sort l0,l1 asc
    if (abs(l0)>abs(l1)) { a=l0; l0=l1; l1=a; }
    if (l0<0.0)          { a=l0; l0=l1; l1=a; }
    if (l0<0.0) return false;
    // relative shift to preserve accuracy after conversion back float
    if (l0>view_depth_max){ a=floor(l0/view_depth_max)*view_depth_max; a-=view_depth_max; view_depth_ll=float(a); if (l1>l0) l1-=a; l0-=a; }
    // conversion back float
    view_depth_l0=float(l0);
    view_depth_l1=float(l1);
    return true;
    }

Portar o restante do shader para dobrar não tem efeito. A única coisa que poderia melhorar isso édouble dados de entrada (entrada édouble mas GL convertê-lo parafloat), mas o GLSL HW atual da mina não permite64 bit interpoladores

Q2 Existe uma maneira de passardouble interpoladores do vértice ao fragment shader?

Algo comovarying dvec4 pixel_pos; no estilo antigo GLSL ouout smooth dvec4 pixel_pos; no perfil principal

questionAnswers(0)

yourAnswerToTheQuestion