cytometer/utils.py
def principal_curvatures_range_image(img, sigma=10):
    """
    Compute Gaussian, Mean and principal curvatures of an image with depth values. Examples of such images
    are topographic maps, range images, depth maps or distance transformations.

    Any of this images can be projected as a Monge patch, a 2D surface embedded in 3D space, f:U->R^3,
    f(x,y) = (x, y, img(x, y)).

    We use a cubic B-spline in tensor-product form to represent the Monge patch by fitting it to the image.
    The reason is that B-splines can be fitted efficiently (scipy.interpolate.RectBivariateSpline), they
    feature compact support, and as they are polynomial functions, the first and second derivatives can be
    easily computed.

    Using these derivatives, there are formulas from elementary Differential Geometry to compute the
    Gaussian, Normal and principal curvatures of the Monge patch.

    Note that these formulas are very sensitive to noise. Thus, we smooth the input image before fitting
    the spline.

    :param img: 2D numpy.ndarray with the distance/depth/range values. The 2D image corresponds to a Monge
    patch, or 2D surface embedded in 3D space of the form (x, y, img(x, y))
    :param sigma: (def sigma=10) Standard deviation in pixels of Gaussian low-pass filtering of the image.
    For sigma=0, no smoothing is performed.
    :return: K, H, k1, k2 = Gaussian curvature, Mean curvature, principal curvature 1, principal
    curvature 2. Each output is an array of the same size as img, with a curvature value per pixel.
    """

    # Note: here x -> rows, y -> columns

    # low-pass filtering
    img = gaussian_filter(img, sigma=sigma)

    # quadratic B-spline interpolation of the image
    sp = RectBivariateSpline(range(img.shape[0]), range(img.shape[1]), img, kx=3, ky=3, s=0)

    # image derivatives and second derivatives
    hx = sp(range(img.shape[0]), range(img.shape[1]), dx=1, grid=True)
    hy = sp(range(img.shape[0]), range(img.shape[1]), dy=1, grid=True)
    hxx = sp(range(img.shape[0]), range(img.shape[1]), dx=2, grid=True)
    hxy = sp(range(img.shape[0]), range(img.shape[1]), dx=1, dy=1, grid=True)
    hyy = sp(range(img.shape[0]), range(img.shape[1]), dy=2, grid=True)

    if DEBUG:
        plt.clf()
        plt.subplot(421)
        plt.imshow(img)
        plt.title('img')
        plt.xlabel('y')
        plt.ylabel('x')
        plt.subplot(422)
        img_resampled = sp(range(img.shape[0]), range(img.shape[1]), grid=True)
        plt.imshow(img_resampled)
        plt.title('B-spline(img)')
        plt.subplot(423)
        plt.imshow(hx)
        plt.title('hx')
        plt.subplot(424)
        plt.imshow(hy)
        plt.title('hy')
        plt.subplot(425)
        plt.imshow(hxx)
        plt.title('hxx')
        plt.subplot(426)
        plt.imshow(hyy)
        plt.title('hyy')
        plt.subplot(427)
        plt.imshow(hxy)
        plt.title('hxy')
        plt.subplot(427)
        plt.imshow(hxy)
        plt.title('hxy')

    # auxiliary products
    hx2 = hx * hx
    hy2 = hy * hy

    # Gaussian curvature
    K = (hxx * hyy - hxy * hxy) / (1 + hx2 + hy2)**2

    # mean curvature
    H = (1 + hx2) * hyy + (1 + hy2) * hxx - 2 * hx * hy * hxy
    H /= 2 * (1 + hx2 + hy2)**1.5

    # principal curvatures
    H2 = H * H
    k1 = H + np.sqrt(H2 - K)
    k2 = H - np.sqrt(H2 - K)

    if DEBUG:
        plt.clf()
        plt.subplot(321)
        plt.imshow(img)
        plt.title('img')
        plt.xlabel('y')
        plt.ylabel('x')
        plt.subplot(323)
        plt.imshow(K)
        plt.title('Gaussian curvature (K)')
        plt.subplot(324)
        plt.imshow(H)
        plt.title('Mean curvature (H)')
        plt.subplot(325)
        plt.imshow(k1.astype(np.float32))
        plt.title('Principal curvature k1')
        plt.subplot(326)
        plt.imshow(k2.astype(np.float32))
        plt.title('Principal curvature k2')

    return K, H, k1, k2

    Compute Gaussian, Mean and principal curvatures of an image with depth values. Examples of such images
    are topographic maps, range images, depth maps or distance transformations.

    Any of this images can be projected as a Monge patch, a 2D surface embedded in 3D space, f:U->R^3,
    f(x,y) = (x, y, img(x, y)).

    We use a cubic B-spline in tensor-product form to represent the Monge patch by fitting it to the image.
    The reason is that B-splines can be fitted efficiently (scipy.interpolate.RectBivariateSpline), they
    feature compact support, and as they are polynomial functions, the first and second derivatives can be
    easily computed.

    Using these derivatives, there are formulas from elementary Differential Geometry to compute the
    Gaussian, Normal and principal curvatures of the Monge patch.

    Note that these formulas are very sensitive to noise. Thus, we smooth the input image before fitting
    the spline.

    :param img: 2D numpy.ndarray with the distance/depth/range values. The 2D image corresponds to a Monge
    patch, or 2D surface embedded in 3D space of the form (x, y, img(x, y))
    :param sigma: (def sigma=10) Standard deviation in pixels of Gaussian low-pass filtering of the image.
    For sigma=0, no smoothing is performed.
    :return: K, H, k1, k2 = Gaussian curvature, Mean curvature, principal curvature 1, principal
    curvature 2. Each output is an array of the same size as img, with a curvature value per pixel.
    
Compute Gaussian, Mean and principal curvatures of an image with depth values. Examples of such images
are topographic maps, range images, depth maps or distance transformations.

Any of this images can be projected as a Monge patch, a 2D surface embedded in 3D space, f:U->R^3,
f(x,y) = (x, y, img(x, y)).

We use a cubic B-spline in tensor-product form to represent the Monge patch by fitting it to the image.
The reason is that B-splines can be fitted efficiently (scipy.interpolate.RectBivariateSpline), they
feature compact support, and as they are polynomial functions, the first and second derivatives can be
easily computed.

Using these derivatives, there are formulas from elementary Differential Geometry to compute the
Gaussian, Normal and principal curvatures of the Monge patch.

Note that these formulas are very sensitive to noise. Thus, we smooth the input image before fitting
the spline.
Compute Gaussian Mean and principal curvatures of an image with depth values . Examples of such images are topographic maps range images depth maps or distance transformations . Any of this images can be projected as a Monge patch a 2D surface embedded in 3D space f : U - > R^3 f ( x y ) = ( x y img ( x y )) . We use a cubic B - spline in tensor - product form to represent the Monge patch by fitting it to the image . The reason is that B - splines can be fitted efficiently ( scipy . interpolate . RectBivariateSpline ) they feature compact support and as they are polynomial functions the first and second derivatives can be easily computed . Using these derivatives there are formulas from elementary Differential Geometry to compute the Gaussian Normal and principal curvatures of the Monge patch . Note that these formulas are very sensitive to noise . Thus we smooth the input image before fitting the spline .
def principal_curvatures_range_image(img, sigma=10):
    img = gaussian_filter(img, sigma=sigma)
    sp = RectBivariateSpline(range(img.shape[0]), range(img.shape[1]), img, kx=3, ky=3, s=0)
    hx = sp(range(img.shape[0]), range(img.shape[1]), dx=1, grid=True)
    hy = sp(range(img.shape[0]), range(img.shape[1]), dy=1, grid=True)
    hxx = sp(range(img.shape[0]), range(img.shape[1]), dx=2, grid=True)
    hxy = sp(range(img.shape[0]), range(img.shape[1]), dx=1, dy=1, grid=True)
    hyy = sp(range(img.shape[0]), range(img.shape[1]), dy=2, grid=True)
    if DEBUG:
        plt.clf()
        plt.subplot(421)
        plt.imshow(img)
        plt.title('img')
        plt.xlabel('y')
        plt.ylabel('x')
        plt.subplot(422)
        img_resampled = sp(range(img.shape[0]), range(img.shape[1]), grid=True)
        plt.imshow(img_resampled)
        plt.title('B-spline(img)')
        plt.subplot(423)
        plt.imshow(hx)
        plt.title('hx')
        plt.subplot(424)
        plt.imshow(hy)
        plt.title('hy')
        plt.subplot(425)
        plt.imshow(hxx)
        plt.title('hxx')
        plt.subplot(426)
        plt.imshow(hyy)
        plt.title('hyy')
        plt.subplot(427)
        plt.imshow(hxy)
        plt.title('hxy')
        plt.subplot(427)
        plt.imshow(hxy)
        plt.title('hxy')
    hx2 = hx * hx
    hy2 = hy * hy
    K = (hxx * hyy - hxy * hxy) / (1 + hx2 + hy2)**2
    H = (1 + hx2) * hyy + (1 + hy2) * hxx - 2 * hx * hy * hxy
    H /= 2 * (1 + hx2 + hy2)**1.5
    H2 = H * H
    k1 = H + np.sqrt(H2 - K)
    k2 = H - np.sqrt(H2 - K)
    if DEBUG:
        plt.clf()
        plt.subplot(321)
        plt.imshow(img)
        plt.title('img')
        plt.xlabel('y')
        plt.ylabel('x')
        plt.subplot(323)
        plt.imshow(K)
        plt.title('Gaussian curvature (K)')
        plt.subplot(324)
        plt.imshow(H)
        plt.title('Mean curvature (H)')
        plt.subplot(325)
        plt.imshow(k1.astype(np.float32))
        plt.title('Principal curvature k1')
        plt.subplot(326)
        plt.imshow(k2.astype(np.float32))
        plt.title('Principal curvature k2')
    return K, H, k1, k2
def principal_curvatures_range_image ( img , sigma = 10 ) : img = gaussian_filter ( img , sigma = sigma ) sp = RectBivariateSpline ( range ( img . shape [ 0 ] ) , range ( img . shape [ 1 ] ) , img , kx = 3 , ky = 3 , s = 0 ) hx = sp ( range ( img . shape [ 0 ] ) , range ( img . shape [ 1 ] ) , dx = 1 , grid = True ) hy = sp ( range ( img . shape [ 0 ] ) , range ( img . shape [ 1 ] ) , dy = 1 , grid = True ) hxx = sp ( range ( img . shape [ 0 ] ) , range ( img . shape [ 1 ] ) , dx = 2 , grid = True ) hxy = sp ( range ( img . shape [ 0 ] ) , range ( img . shape [ 1 ] ) , dx = 1 , dy = 1 , grid = True ) hyy = sp ( range ( img . shape [ 0 ] ) , range ( img . shape [ 1 ] ) , dy = 2 , grid = True ) if DEBUG : plt . clf ( ) plt . subplot ( 421 ) plt . imshow ( img ) plt . title ( 'img' ) plt . xlabel ( 'y' ) plt . ylabel ( 'x' ) plt . subplot ( 422 ) img_resampled = sp ( range ( img . shape [ 0 ] ) , range ( img . shape [ 1 ] ) , grid = True ) plt . imshow ( img_resampled ) plt . title ( 'B-spline(img)' ) plt . subplot ( 423 ) plt . imshow ( hx ) plt . title ( 'hx' ) plt . subplot ( 424 ) plt . imshow ( hy ) plt . title ( 'hy' ) plt . subplot ( 425 ) plt . imshow ( hxx ) plt . title ( 'hxx' ) plt . subplot ( 426 ) plt . imshow ( hyy ) plt . title ( 'hyy' ) plt . subplot ( 427 ) plt . imshow ( hxy ) plt . title ( 'hxy' ) plt . subplot ( 427 ) plt . imshow ( hxy ) plt . title ( 'hxy' ) hx2 = hx * hx hy2 = hy * hy K = ( hxx * hyy - hxy * hxy ) / ( 1 + hx2 + hy2 ) ** 2 H = ( 1 + hx2 ) * hyy + ( 1 + hy2 ) * hxx - 2 * hx * hy * hxy H /= 2 * ( 1 + hx2 + hy2 ) ** 1.5 H2 = H * H k1 = H + np . sqrt ( H2 - K ) k2 = H - np . sqrt ( H2 - K ) if DEBUG : plt . clf ( ) plt . subplot ( 321 ) plt . imshow ( img ) plt . title ( 'img' ) plt . xlabel ( 'y' ) plt . ylabel ( 'x' ) plt . subplot ( 323 ) plt . imshow ( K ) plt . title ( 'Gaussian curvature (K)' ) plt . subplot ( 324 ) plt . imshow ( H ) plt . title ( 'Mean curvature (H)' ) plt . subplot ( 325 ) plt . imshow ( k1 . astype ( np . float32 ) ) plt . title ( 'Principal curvature k1' ) plt . subplot ( 326 ) plt . imshow ( k2 . astype ( np . float32 ) ) plt . title ( 'Principal curvature k2' ) return K , H , k1 , k2
Compute Gaussian, Mean and principal curvatures of an image with depth values.
Compute Gaussian Mean and principal curvatures of an image with depth values .
"""
    Compute Gaussian, Mean and principal curvatures of an image with depth values. Examples of such images
    are topographic maps, range images, depth maps or distance transformations.

    Any of this images can be projected as a Monge patch, a 2D surface embedded in 3D space, f:U->R^3,
    f(x,y) = (x, y, img(x, y)).

    We use a cubic B-spline in tensor-product form to represent the Monge patch by fitting it to the image.
    The reason is that B-splines can be fitted efficiently (scipy.interpolate.RectBivariateSpline), they
    feature compact support, and as they are polynomial functions, the first and second derivatives can be
    easily computed.

    Using these derivatives, there are formulas from elementary Differential Geometry to compute the
    Gaussian, Normal and principal curvatures of the Monge patch.

    Note that these formulas are very sensitive to noise. Thus, we smooth the input image before fitting
    the spline.

    :param img: 2D numpy.ndarray with the distance/depth/range values. The 2D image corresponds to a Monge
    patch, or 2D surface embedded in 3D space of the form (x, y, img(x, y))
    :param sigma: (def sigma=10) Standard deviation in pixels of Gaussian low-pass filtering of the image.
    For sigma=0, no smoothing is performed.
    :return: K, H, k1, k2 = Gaussian curvature, Mean curvature, principal curvature 1, principal
    curvature 2. Each output is an array of the same size as img, with a curvature value per pixel.
    """ # Note: here x -> rows, y -> columns # low-pass filtering # quadratic B-spline interpolation of the image # image derivatives and second derivatives # auxiliary products # Gaussian curvature # mean curvature # principal curvatures
K, H, k1, k2 = Gaussian curvature, Mean curvature, principal curvature 1, principal
curvature 2. Each output is an array of the same size as img, with a curvature value per pixel.
K H k1 k2 = Gaussian curvature Mean curvature principal curvature 1 principal curvature 2 . Each output is an array of the same size as img with a curvature value per pixel .
img
2D numpy.ndarray with the distance/depth/range values. The 2D image corresponds to a Monge
patch, or 2D surface embedded in 3D space of the form (x, y, img(x, y))
2D numpy . ndarray with the distance / depth / range values . The 2D image corresponds to a Monge patch or 2D surface embedded in 3D space of the form ( x y img ( x y ))
sigma
(def sigma=10) Standard deviation in pixels of Gaussian low-pass filtering of the image.
For sigma=0, no smoothing is performed.
( def sigma = 10 ) Standard deviation in pixels of Gaussian low - pass filtering of the image . For sigma = 0 no smoothing is performed .