Минус значения коэффициентов искажения приводят к искажению изображений с помощью OpenGL Vertex Shader

Question

Я работаю с python moderngl для генерации имитированных видов сцены из облаков точек, сделанных из цифровой модели рельефа и бортовых фотографий. Я реализовал модель камеры на основе OpenCV в GLSL Vertex Shader, и она обычно работает хорошо.

Однако некоторые минусовые значения коэффициентов радиального искажения (k1 или k2) приводят к искажению изображений.

Здесь ниже это код вершинного шейдера и фрагментного шейдера.

    prog = ctx.program(
     vertex_shader='''
         #version 330
         in vec3 in_vert;
         in vec3 in_color;
         out vec4 v_color;

         // decrare some values used inside GPU by "uniform"
         // the real values will be later set by CPU
         uniform mat4 proj; // projection matrix
         uniform mat4 view; // model view matrix
         uniform mat4 rotate; // rotation matrix for normals
         uniform float dist_coeffs[5]; // distortion coefficients 
         // distortion coefficients. k1x, k1y, k2x, k2y, p1, p2
         // radial distortion parameters are divided as x-axis and y-axis
         vec4 distort(vec4 view_pos){
          // normalize
          float z =  view_pos[2];
          float z_inv = 1 / z;
          float x1 = view_pos[0] * z_inv;
          float y1 = view_pos[1] * z_inv;

          // precalculations
          float x1_2 = x1*x1;
          float y1_2 = y1*y1;
          float x1_y1 = x1*y1;
          float r2 = x1_2 + y1_2;
          float r4 = r2*r2;

          // radial distortion factor
          float r_dist_x = (1 + dist_coeffs[0] * r2 + dist_coeffs[2] * r4);
          float r_dist_y = (1 + dist_coeffs[1] * r2 + dist_coeffs[3] * r4);

          // full (radial + tangential) distortion
          float x2 = x1*r_dist_x + 2*dist_coeffs[4]*x1_y1 + dist_coeffs[5]*(r2 + 2*x1_2);
          float y2 = y1*r_dist_y + 2*dist_coeffs[5]*x1_y1 + dist_coeffs[4]*(r2 + 2*y1_2);

          // denormalize for projection
          return vec4(x2*z, y2*z, z, view_pos[3]);
          }

         void main() {
             vec4 local_pos = vec4(in_vert, 1.0);
             vec4 view_pos = vec4(view * local_pos);
             vec4 dist_pos = distort(view_pos);
             vec4 v_norm = rotate * vec4(in_color, 1.0); // normals are transformed to camera CRS, -1~1
             v_color = (v_norm + 1) / 2; // normals as rgb, 0~1
             gl_Position = vec4(proj * dist_pos);
         }
     ''',
     fragment_shader='''
         #version 330
         in vec4 v_color;
         layout(location=0)out vec4 f_color;
         void main() {
             f_color = vec4(v_color); // 1,0 added is alpha
         }
     '''
    )

А здесь ниже я показываю программу R, которая создает вид модели и матрицу перспективной проекции.

# camera intrinsic matrix
projection_mat <- function(fov_x_degree, w, h, near = -1, far = 1
                           , cx = w/2, cy = h/2){
  fov_x <- fov_x_degree * pi /180
  fov_y <- fov_x * h / w
  fx <- 1/(tan(fov_x/2))
  fy <- 1/(tan(fov_y/2))
  matrix(c(fx , 0, (w-2*cx)/w, 0, 
           0, fy, -(h-2*cy)/h, 0,
           0, 0, -(far+near)/(far-near), -2*far*near/(far-near),
           0, 0, -1, 0), 
         ncol = 4, 
         byrow = F) %>% 
    as.vector() %>% 
    reticulate::np_array(dtype = "float32")
}

# camera extrinsic parameter
modelview_mat <- function(pan_degree, tilt_degree, roll_degree, t_x, t_y, t_z){
  y <-  (360 - pan_degree) * pi / 180
  x <-  (tilt_degree) * pi / 180 
  z <-  (roll_degree) * pi / 180

  rmat_z <- c(cos(z), -sin(z), 0, 0,
              sin(z), cos(z), 0, 0,
              0, 0, 1, 0,
              0, 0, 0, 1) %>% matrix(ncol = 4, byrow = T)
  rmat_x <- c(1, 0, 0, 0,
              0, cos(x), -sin(x), 0,
              0, sin(x), cos(x), 0,
              0, 0, 0, 1) %>% matrix(ncol = 4, byrow = T)
  rmat_y <- c(cos(y), 0, sin(y), 0,
              0, 1, 0, 0,
              -sin(y), 0, cos(y), 0,
              0, 0, 0, 1) %>%matrix(ncol = 4, byrow = T)
  rmat <- rmat_x %*% rmat_z %*% rmat_y
  tmat <- matrix(c(1, 0, 0, -t_x,
                   0, 1, 0, -t_z,
                   0, 0, 1, -t_y,
                   0, 0, 0, 1), byrow = T, ncol = 4)
  rmat %*% tmat %>%
    as.vector() %>% 
    reticulate::np_array(dtype = "float32")
}

это сгенерированное изображение без ошибок

и это сломанный