Looking for an OpenCV4 alternative for OpenCV2 Interpolation Flag/Remapping Function

Question

I am currently porting a Python 2.7 code written around 10 years ago to Python 3, and it was written using OpenCV 2.x and in one line there is a part that warps an image, however when I try to run that line in Python 3 using OpenCV 4.x, I cannot run it due to the interpolation flags not being recognized by OpenCV.

Here is the line:

cv.Remap(im, output_image, ix, iy, flags=cv.CV_INTER_LINEAR+cv.CV_WARP_FILL_OUTLIERS+cv.CV_WARP_INVERSE_MAP)

There is a very similar function of OpenCV in its current version, and this line without the other interpolation flags work, however the result seems to be different:

output_image = cv.remap(im,output_image,ix,iy,flags=cv.INTER_LINEAR)

im is the input image, output_image is an empty array that has been declared before, and ix and iy are two arrays. The code is a Lightfield Microscopy calibration code, and in this function they are doing a "cubic warp", which I couldn't also find any information on the internet.

I tried using warpPerspective, but it is asking for a transformation matrix, which I don't have. I tried using findTransformECC, but couldn't manage to make it work. I tried using OpenCV 2.4 in Python 2.7 and it works, but as Python 2 is not supported anymore, it would be better to have a solution that is more future-proof.

Here is the full function:

# Warp an image.
    #
    # This function expects im to be a numpy float32 array.  Returns 
    def warp_image(self, input_image, output_pixels_per_lenslet, direction="R",
                   cropToInside = False, lenslet_offset = None, output_size = None):
        #im = cv.fromarray(input_image)
        im = (input_image)
        ul = self.eval_point([0, 0], 'f')
        #ur = self.eval_point([0, im.cols], 'f')
        ur = self.eval_point([0, im.shape[1]], 'f')
        #ll = self.eval_point([im.rows, 0], 'f')
        ll = self.eval_point([im.shape[0],0], 'f')
        #lr = self.eval_point([im.rows, im.cols], 'f')
        lr = self.eval_point([im.shape[0], im.shape[1]], 'f')
        
        leftbound = np.ceil(max(ul[1], ll[1]))
        rightbound = np.floor(min(ur[1], lr[1]))
        topbound = np.ceil(max(ul[0], ur[0]))
        bottombound = np.floor(min(ll[0], lr[0]))

        # Don't crop left of lenslets (0, y) or above (x, 0)
        leftbound = max(leftbound, 0)
        topbound = max(topbound,0)

        if output_size != None:
            putative_output_size = output_size
        else:
            nt = int(np.floor(bottombound - topbound))
            ns = int(np.floor(rightbound - leftbound))
            putative_output_size = (nt*output_pixels_per_lenslet, ns*output_pixels_per_lenslet)
        
        # Create the output image
        #output_image = cv.CreateMat(putative_output_size[0], putative_output_size[1], im.type)
        output_image = np.empty(shape=(putative_output_size[0],putative_output_size[1]),dtype = im.dtype)
        # Apply the transform.
        scaled_shift = (0.0, 0.0)
        if (direction == 'f' or direction == 'F'):
            coeff = np.copy(self.forwardCoefficients)
            coeff[:,1:] *= output_pixels_per_lenslet
            coeff[:,3:] *= output_pixels_per_lenslet
            coeff[:,6:] *= output_pixels_per_lenslet
            if lenslet_offset != None:
                scaled_shift = (lenslet_offset[0] / output_pixels_per_lenslet,
                                lenslet_offset[1] / output_pixels_per_lenslet)
        else:
            coeff = np.copy(self.reverseCoefficients)
            coeff[:,1:] /= output_pixels_per_lenslet
            coeff[:,3:] /= output_pixels_per_lenslet
            coeff[:,6:] /= output_pixels_per_lenslet
            if np.any(lenslet_offset) != None:
                scaled_shift = (lenslet_offset[0] * output_pixels_per_lenslet,
                                lenslet_offset[1] * output_pixels_per_lenslet)

        (x, y) = np.meshgrid(np.arange(putative_output_size[1], dtype='float32'), np.arange(putative_output_size[0], dtype='float32'))

        ix_array = (np.ones((putative_output_size[0], putative_output_size[1])).astype('float32') * coeff[0,0] + x * coeff[0,1] + y * coeff[0,2] + x * x * coeff[0,3] + x * y * coeff[0,4] + y * y * coeff[0,5] + x * x * x * coeff[0,6] + x * x * y * coeff[0,7] + x * y * y * coeff[0,8] + y * y * y * coeff[0,9] + scaled_shift[0])
        iy_array = (np.ones((putative_output_size[0], putative_output_size[1])).astype('float32') * coeff[1,0] + x * coeff[1,1] + y * coeff[1,2] + x * x * coeff[1,3] + x * y * coeff[1,4] + y * y * coeff[1,5] + x * x * x * coeff[1,6] + x * x * y * coeff[1,7] + x * y * y * coeff[1,8] + y * y * y * coeff[1,9] + scaled_shift[1])
        #ix = cv.fromarray(ix_array.astype(np.float32))
        #iy = cv.fromarray(iy_array.astype(np.float32))
        ix = ix_array.astype(np.float32)
        iy = iy_array.astype(np.float32)

        
        cv.Remap(im, output_image, ix, iy, flags=cv.CV_INTER_LINEAR+cv.CV_WARP_FILL_OUTLIERS+cv.CV_WARP_INVERSE_MAP)
        result = np.asarray(output_image)

        # Don't crop left of lenslets (0, y) or above (x, 0)
        leftbound = max(leftbound, 0)
        topbound = max(topbound,0)
    
        if cropToInside:
            return result[int(topbound * output_pixels_per_lenslet):int(bottombound * output_pixels_per_lenslet),
                          int(leftbound * output_pixels_per_lenslet):int(rightbound * output_pixels_per_lenslet)]
        
        else:
            return result