Otras fuentes parecen dar diferentes definiciones de los núcleos más grandes. La biblioteca Intel IPP, por ejemplo, da el núcleo 5x5 como

1  2 0  -2 -1
4  8 0  -8 -4
6 12 0 -12 -6
4  8 0  -8 -4
1  2 0  -2 -1

Intuitivamente, esto tiene más sentido para mí porque estás prestando más atención a los elementos más cercanos al centro. También tiene una definición natural en términos del núcleo 3x3 que es fácil de extender para generar núcleos más grandes. Dicho esto, en mi breve búsqueda he encontrado 3 definiciones diferentes del núcleo 5x5, así que sospecho que (como Paul dice) los núcleos más grandes son ad hoc, por lo que esta no es de ninguna manera la respuesta definitiva.

El núcleo 3x3 es el producto externo de un núcleo de suavizado y un núcleo de gradiente, en Matlab esto es algo así como

sob3x3 = [ 1 2 1 ]' * [1 0 -1]

Los núcleos más grandes se pueden definir convolviendo el núcleo 3x3 con otro núcleo de suavizado

sob5x5 = conv2( [ 1 2 1 ]' * [1 2 1], sob3x3 )

Puede repetir el proceso para obtener núcleos progresivamente más grandes

sob7x7 = conv2( [ 1 2 1 ]' * [1 2 1], sob5x5 )
sob9x9 = conv2( [ 1 2 1 ]' * [1 2 1], sob7x7 )
...

Hay muchas otras formas de escribirlo, pero creo que esto explica exactamente lo que está sucediendo mejor. Básicamente, comienzas con un núcleo suavizado en una dirección y una estimación de diferencias finitas de la derivada en la otra y luego solo aplicas suavizado hasta obtener el tamaño del núcleo que deseas.

Debido a que es solo una serie de convoluciones, todas las propiedades agradables tienen, (conmutatividad, asociatividad, etc.) que podrían ser útiles para su implementación. Por ejemplo, puede separar trivialmente el núcleo 5x5 en su suavizado y derivado componentes:

Sob5x5 = conv([1 2 1],[1 2 1])' * conv([1 2 1],[-1 0 1])

Tenga en cuenta que para ser un estimador derivado "adecuado", el Sobel 3x3 debe escalarse por un factor de 1/8:

sob3x3 = 1/8 * [ 1 2 1 ]' * [1 0 -1]

Y cada núcleo más grande necesita ser escalado por un factor adicional de 1/16 (porque los núcleos suavizados no están normalizados):

sob5x5 = 1/16 * conv2( [ 1 2 1 ]' * [1 2 1], sob3x3 )
sob7x7 = 1/16 * conv2( [ 1 2 1 ]' * [1 2 1], sob5x5 )
...

Solución completa para tamaños y ángulos arbitrarios del kernel de Sobel

Tl; dr: salta a la sección'Ejemplos'

Para agregar otra solución, expandiendo este documento (no es particularmente de alta calidad, pero muestra algunos gráficos y matrices utilizables que comienzan en la parte inferior de la página 2).

Objetivo

Lo que estamos tratando de hacer es estimar el gradiente local de la imagen en la posición (x,y). El gradiente es un vector formado por los componentes en dirección x e y, gx y gy.

Ahora, imagine que queremos aproximar el gradiente basado en nuestro píxel (x,y) y sus vecinos como una operación del núcleo (3x3, 5x5, o cualquier tamaño).

Idea de solución

Podemos aproximar el gradiente sumando las proyecciones de todos los pares vecino-centro en la dirección del gradiente. (El núcleo de Sobel es solo un método particular de ponderación de las diferentes contribuciones, y también lo es Prewitt, básicamente).

Explícito pasos intermedios para 3x3

Esta es la imagen local, píxel central (x,y) marcado como' o ' (centro)

a b c
d o f
g h i

Digamos que queremos el gradiente en dirección x positiva. El vector unitario en la dirección x positiva es (1,0) [Luego usaré la convención de que la dirección y positiva es HACIA ABAJO, es decir, (0,1), y que (0,0) es la parte superior izquierda de la imagen).]

El vector de o a f ('of' para abreviar) es (1,0). El gradiente en la dirección ' de ' es (f-o) / 1 (valor de la imagen en el píxel aquí denotado f menos el valor en el centro o, dividido por la distancia entre esos píxeles). Si proyectamos el vector unitario de ese gradiente vecino particular en nuestra dirección de gradiente deseada (1,0) a través de un producto escalar obtenemos 1. Aquí hay una pequeña tabla con las contribuciones de todos los vecinos, comenzando con los casos más fáciles. Tenga en cuenta que para las diagonales, su distancia es sqrt2, y los vectores unitarios en las direcciones diagonales son 1 / sqrt2 * (+/-1, +/-1)

f:   (f-o)/1     * 1
d:   (d-o)/1     * -1       because (-1, 0) dot (1, 0) = -1
b:   (b-o)/1     * 0        because (0, -1) dot (1, 0) = 0
h:   (h-o)/1     * 0        (as per b)
a:   (a-o)/sqrt2 * -1/sqrt2 distance is sqrt2, and 1/sqrt2*(-1,-1) dot (1,0) = -1/sqrt2
c:   (c-o)/sqrt2 * +1/sqrt2   ...
g:   (g-o)/sqrt2 * -1/sqrt2   ...
i:   (i-o)/sqrt2 * +1/sqrt2   ...

Editar para clarificar: Hay dos factores de 1 / sqrt (2) por la siguiente razón:

1. Estamos interesados en la contribución al gradiente en una dirección específica (aquí x), por lo que necesitamos proyectar el gradiente direccional desde el píxel central al píxel vecino en la dirección que nos interesa. Esto se logra tomando el producto escalar de los vectores unitarios en las direcciones respectivas, lo que introduce el primer factor 1 / L(aquí 1/sqrt (2) para el diagonales).

2. El gradiente mide el cambio infinitesimal en un punto, que aproximamos por diferencias finitas. En términos de una ecuación lineal, m = (y2-y1)/(x2-x1). Por esta razón, la diferencia de valor entre el píxel central y el píxel vecino (y2-y1) debe distribuirse sobre su distancia (corresponde a x2-x1) para obtener las unidades de ascenso por unidad de distancia. Esto produce un segundo factor de 1 / L(aquí 1/sqrt (2) para las diagonales)

Ok, ahora sabemos las contribuciones. Simplifiquemos esta expresión combinando pares opuestos de contribuciones de píxeles. Voy a empezar con d y f:

{(f-o)/1 * 1} + {(d-o)/1 * -1}
= f - o - (d - o)
= f - d

Ahora la primera diagonal:

{(c-o)/sqrt2 * 1/sqrt2} + {(g-o)/sqrt2 * -1/sqrt2}
= (c - o)/2 - (g - o)/2
= (c - g)/2

La segunda diagonal contribuye (i - a)/2. La dirección perpendicular contribuye a cero. Ten en cuenta que todas las contribuciones del píxel central 'o' desaparecen.

Ahora hemos calculado las contribuciones de todos los vecinos más cercanos al gradiente en dirección x positiva en pixel (x, y), por lo que nuestro la aproximación total del gradiente en la dirección x es simplemente su suma:

gx(x,y) = f - d + (c - g)/2 + (i - a)/2

Podemos obtener el mismo resultado usando un núcleo de convolución donde los coeficientes se escriben en el lugar del píxel vecino correspondiente:

-1/2  0  1/2
 -1   0   1
-1/2  0  1/2

Si no quiere tratar con fracciones, multiplique esto por 2 y obtenga el conocido núcleo Sobel 3x3.

      -1 0 1
G_x = -2 0 2
      -1 0 1

La multiplicación por dos solo sirve para obtener números enteros convenientes. El escalado de la imagen de salida es básicamente arbitrario, la mayoría de las veces lo normaliza a su rango de imagen, de todos modos (para obtener resultados claramente visibles).

Por el mismo razonamiento que el anterior, se obtiene el núcleo para el gradiente vertical gy proyectando las contribuciones vecinas sobre el vector unitario en dirección y positiva (0,1)

      -1 -2 -1
G_y =  0  0  0
       1  2  1

Fórmula para núcleos de tamaño arbitrario

Si desea núcleos 5x5 o más grandes, solo necesita prestar atención a las distancias, por ejemplo,

A B 2 B A
B C 1 C B
2 1 - 1 2
B C 1 C B
A B 2 B A

Donde

A = 2 * sqrt2
B = sqrt5
C = sqrt2.

Si la longitud del vector que conecta dos píxeles cualesquiera es L, el vector unitario en esa dirección tiene un prefactor de 1 / L. Por esta razón, las contribuciones de cualquier píxel ' k 'a (digamos) el gradiente x (1,0) se pueden simplificar a" (diferencia de valor sobre la distancia cuadrada) veces (DotProduct de vector de dirección no normalizado ' ok ' con vector de gradiente, por ejemplo (1,0) )"

gx_k = (k - o)/(pixel distance^2) ['ok' dot (1,0)].

Debido a que el producto escalar del vector de conexión con el vector unitario x selecciona el entrada vectorial, la entrada correspondiente del kernel G_x en la posición k es simplemente

i / (i*i + j*j)

Donde i y j son el número de pasos desde el píxel central hasta el píxel k en dirección x e y. En el cálculo 3x3 anterior, el píxel ' a 'tendría i = -1 (1 a la izquierda), j = -1 (1 a la parte superior) y, por lo tanto, la entrada del núcleo' a ' es -1 / (1 + 1) = -1/2.

Las entradas para el núcleo G_y son

j/(i*i + j*j). 

Si quiero valores enteros para mi núcleo, sigo estos pasos:

• compruebe el rango disponible de la imagen de salida
• calcule el resultado más alto posible de aplicar el kernel de coma flotante (es decir, asuma el valor máximo de entrada bajo todas las entradas positivas del kernel, por lo que el valor de salida es (suma sobre todos los valores positivos del kernel) * (valor máximo de la imagen de entrada posible). Si ha firmado la entrada, también debe considerar los valores negativos. El peor caso es entonces la suma de todos los valores positivos + suma de todos los valores abs de entradas negativas (si entrada máxima en positivos, - entrada máxima en negativos). editar: la suma de todos los valores de abs también se ha llamado acertadamente el peso del núcleo
• calcular el escalado máximo permitido para el núcleo (sin rebosar el rango de la imagen de salida)
• para todos los múltiplos enteros (desde 2 hasta un máximo superior) del núcleo de coma flotante: compruebe cuál tiene la suma más baja de errores de redondeo absolutos y use este núcleo

Así que en resumen:

Gx_ij = i / (i*i + j*j)
Gy_ij = j / (i*i + j*j)

Donde i, j es posición en el núcleo contada desde el centro. Escala las entradas del núcleo según sea necesario para obtener números enteros (o al menos aproximaciones cercanas).

Estas fórmulas se aplican a todos los tamaños de kernel.

Ejemplos

          -2/8 -1/5  0  1/5  2/8           -5  -4  0   4   5
          -2/5 -1/2  0  1/2  2/5           -8 -10  0  10   8
G_x (5x5) -2/4 -1/1  0  1/1  2/4  (*20) = -10 -20  0  20  10
          -2/5 -1/2  0  1/2  2/5           -8 -10  0  10   8
          -2/8 -1/5  0  1/5  2/8           -5  -4  0   4   5

Tenga en cuenta que los píxeles centrales de 3x3 del núcleo de 5x5 en notación flotante son solo el núcleo de 3x3, es decir, los núcleos más grandes representan una aproximación continua con datos adicionales pero de menor peso. Esto continúa en el núcleo más grande tamaños:

           -3/18 -2/13 -1/10 0  1/10 2/13 3/18
           -3/13 -2/8  -1/5  0  1/5  2/8  3/13
           -3/10 -2/5  -1/2  0  1/2  2/5  3/10
G_x (7x7)  -3/9  -2/4  -1/1  0  1/1  2/4  3/9 
           -3/10 -2/5  -1/2  0  1/2  2/5  3/10
           -3/13 -2/8  -1/5  0  1/5  2/8  3/13
           -3/18 -2/13 -1/10 0  1/10 2/13 3/18

Las representaciones exactas de enteros se vuelven imprácticas en este punto.

Por lo que puedo decir (no tengo acceso al documento original), la parte de "Sobel" de esto es ponderar adecuadamente las contribuciones. La solución Prewitt se puede obtener omitiendo la ponderación de distancia e introduciendo i y j en el núcleo según corresponda.

Bonus: Núcleos Sobel para direcciones arbitrarias

Para que podamos aproximar los componentes x e y de la imagen gradiente (que en realidad es un vector, como se indica al principio). El gradiente en cualquier dirección arbitraria alfa (medido matemáticamente positivo, en este caso en el sentido de las agujas del reloj ya que positivo y es hacia abajo) se puede obtener proyectando el vector de gradiente sobre el vector unitario de gradiente alfa.

El vector alfa-unidad es (cos alfa, sin alfa). Para alfa = 0° se puede obtener el resultado para gx, para alfa = 90° se obtiene gy.

g_alpha = (alpha-unit vector) dot (gx, gy)
        = (cos a, sin a) dot (gx, gy)
        = cos a * gx + sin a * gy

Si se molesta en escribir gx y gy como sumas de las contribuciones vecinas, se da cuenta de que puede agrupar la expresión larga resultante por términos que se aplican al mismo píxel vecino, y luego reescribir esto como un solo núcleo de convolución con entradas

G_alpha_ij = (i * cos a + j * sin a)/(i*i + j*j)

Si desea la aproximación de enteros más cercana, siga los pasos descritos anteriormente.

Rápidamente hackeé un algoritmo para generar un núcleo Sobel de cualquier tamaño impar > 1, basado en los ejemplos dados por @Paul R:

    public static void CreateSobelKernel(int n, ref float[][] Kx, ref float[][] Ky)
    {
        int side = n * 2 + 3;
        int halfSide = side / 2;
        for (int i = 0; i < side; i++)
        {
            int k = (i <= halfSide) ? (halfSide + i) : (side + halfSide - i - 1);
            for (int j = 0; j < side; j++)
            {
                if (j < halfSide)
                    Kx[i][j] = Ky[j][i] = j - k;
                else if (j > halfSide)
                    Kx[i][j] = Ky[j][i] = k - (side - j - 1);
                else
                    Kx[i][j] = Ky[j][i] = 0;
            }
        }
    }

Espero que ayude.

Generador de filtro de gradiente Sobel

(Esta respuesta se refiere al análisis dado por @Daniel, arriba.)

Gx[i,j] = i / (i*i + j*j)

Gy[i,j] = j / (i*i + j*j)

Este es un resultado importante, y una explicación mejor que la que se puede encontrar en el documento original. Debería estar escrito en Wikipedia , o en alguna parte, porque también parece superior a cualquier otra discusión del tema disponible en Internet.

Sin embargo, en realidad no es cierto que las representaciones con valores enteros sean poco práctico para filtros de tamaño superior a 5 * 5, como se afirma. Usando enteros de 64 bits, los tamaños de filtro Sobel de hasta 15 * 15 se pueden expresar exactamente.

Aquí están los cuatro primeros; el resultado debe dividirse por el "peso", de modo que el gradiente de una región de imagen como la siguiente, se normalice a un valor de 1.

1 2 3 4 5
1 2 3 4 5
1 2 3 4 5
1 2 3 4 5
1 2 3 4 5

Gx (3):

-1/2  0/1  1/2           -1  0  1
-1/1    0  1/1   * 2 =   -2  0  2
-1/2  0/1  1/2           -1  0  1

weight = 4               weight = 8

Gx (5):

-2/8 -1/5  0/4  1/5  2/8             -5  -4   0   4   5
-2/5 -1/2  0/1  1/2  2/5             -8 -10   0  10   8
-2/4 -1/1    0  1/1  2/4   * 20 =   -10 -20   0  20  10
-2/5 -1/2  0/1  1/2  2/5             -8 -10   0  10   8
-2/8 -1/5  0/4  1/5  2/8             -5  -4   0   4   5

weight = 12                          weight = 240

Gx (7) :

-3/18 -2/13 -1/10   0/9  1/10  2/13  3/18             -130 -120  -78    0   78  120  130
-3/13  -2/8  -1/5   0/4   1/5   2/8  3/13             -180 -195 -156    0  156  195  180
-3/10  -2/5  -1/2   0/1   1/2   2/5  3/10             -234 -312 -390    0  390  312  234
 -3/9  -2/4  -1/1     0   1/1   2/4   3/9   * 780 =   -260 -390 -780    0  780  390  260
-3/10  -2/5  -1/2   0/1   1/2   2/5  3/10             -234 -312 -390    0  390  312  234
-3/13  -2/8  -1/5   0/4   1/5   2/8  3/13             -180 -195 -156    0  156  195  180
-3/18 -2/13 -1/10   0/9  1/10  2/13  3/18             -130 -120  -78    0   78  120  130

weight = 24                                           weight = 18720

Gx (9) :

-4/32 -3/25 -2/20 -1/17  0/16  1/17  2/20  3/25  4/32                -16575  -15912  -13260   -7800       0    7800   13260   15912   16575
-4/25 -3/18 -2/13 -1/10   0/9  1/10  2/13  3/18  4/25                -21216  -22100  -20400  -13260       0   13260   20400   22100   21216
-4/20 -3/13  -2/8  -1/5   0/4   1/5   2/8  3/13  4/20                -26520  -30600  -33150  -26520       0   26520   33150   30600   26520
-4/17 -3/10  -2/5  -1/2   0/1   1/2   2/5  3/10  4/17                -31200  -39780  -53040  -66300       0   66300   53040   39780   31200
-4/16  -3/9  -2/4  -1/1     0   1/1   2/4   3/9  4/16   * 132600 =   -33150  -44200  -66300 -132600       0  132600   66300   44200   33150
-4/17 -3/10  -2/5  -1/2   0/1   1/2   2/5  3/10  4/17                -31200  -39780  -53040  -66300       0   66300   53040   39780   31200
-4/20 -3/13  -2/8  -1/5   0/4   1/5   2/8  3/13  4/20                -26520  -30600  -33150  -26520       0   26520   33150   30600   26520
-4/25 -3/18 -2/13 -1/10   0/9  1/10  2/13  3/18  4/25                -21216  -22100  -20400  -13260       0   13260   20400   22100   21216
-4/32 -3/25 -2/20 -1/17  0/16  1/17  2/20  3/25  4/32                -16575  -15912  -13260   -7800       0    7800   13260   15912   16575

weight = 40                                                          weight = 5304000

El programa Ruby adjunto a continuación, calculará Los filtros Sobel y los pesos correspondientes de cualquier tamaño, aunque los filtros con valores enteros no son probablemente útiles para tamaños mayores de 15*15.

#!/usr/bin/ruby

# Sobel image gradient filter generator
# by <[email protected]> -- Sept 2017
# reference:
# https://stackoverflow.com/questions/9567882/sobel-filter-kernel-of-large-size


if (s = ARGV[0].to_i) < 3 || (s % 2) == 0
    $stderr.puts "invalid size"
    exit false
end

s /= 2


n = 1

# find least-common-multiple of all fractional denominators
(0..s).each { |j|
    (1..s).each { |i|
        d = i*i + j*j
        n = n.lcm(d / d.gcd(i))
    }
}


fw1 = format("%d/%d", s, 2*s*s).size + 2
fw2 = format("%d", n).size + 2


weight = 0
s1 = ""
s2 = ""

(-s..s).each { |y|
    (-s..s).each { |x|
        i, j = x, y   # "i, j = y, x" for transpose
        d = i*i + j*j
        if (i != 0)
            if (n * i % d) != 0   # this should never happen
                $stderr.puts "inexact division: #{n} * #{i} / ((#{i})^2 + (#{j})^2)"
                exit false
            end
            w = n * i / d
            weight += i * w
        else
            w = 0
        end
        s1 += "%*s" % [fw1, d > 0 ? "%d/%d" % [i, d] : "0"]
        s2 += "%*d" % [fw2, w]
    }
    s1 += "\n" ; s2 += "\n"
}


f = n.gcd(weight)

puts s1

puts "\nweight = %d%s" % [weight/f, f < n ? "/%d" % (n/f) : ""]

puts "\n* #{n} =\n\n"

puts s2

puts "\nweight = #{weight}"

Como Adam Bowen explicó en su respuesta , el núcleo Sobel es una combinación de un suavizado a lo largo de un eje, y una derivada de diferencia central a lo largo del otro eje:{[18]]}

sob3x3 = [1 2 1]' * [1 0 -1]

El suavizado añade regularización (reduce la sensibilidad al ruido).

(Estoy omitiendo todos los factores 1/8 en este post, al igual que el propio Sobel, lo que significa que el operador determina la derivada hasta la escala. También, * siempre significa convolución en este post.)

Generalicemos esto: {[18]]}

deriv_kernel = smoothing_kernel * d/dx

Una de las propiedades de la convolución es que

d/dx f = d/dx * f

Es decir, enrollar una imagen con el operador derivado elemental produce la derivada de la imagen. Observando también que la convolución es conmutativa,

deriv_kernel = d/dx * smoothing_kernel = d/dx smoothing_kernel

Es decir, el núcleo derivado es el derivado de un núcleo de suavizado.

Tenga en cuenta que la aplicación de un núcleo de este tipo a una imagen por convolución:

image * deriv_kernel = image * smoothing_kernel * d/dx = d/dx (image * smoothing_kernel)

Es decir, con esto núcleo derivado idealizado y generalizado podemos calcular el verdadero derivado de la imagen suavizada. Este no es, por supuesto, el caso con el núcleo Sobel, ya que utiliza una aproximación de diferencia central a la derivada. Pero elegir un mejor smoothing_kernel, esto se puede lograr. El núcleo gaussiano es la opción ideal aquí, ya que ofrece el mejor compromiso entre compacidad en el dominio espacial (pequeña huella del núcleo) con compacidad en el dominio de frecuencia (buen suavizado). Además, la gaussiana es perfectamente isotrópica y separable. El uso de un núcleo derivado gaussiano produce el mejor operador derivado regularizado posible.

Por lo tanto, si está buscando un operador Sobel más grande, porque necesita más regularización, use un operador derivado gaussiano en su lugar.

Analicemos un poco más el núcleo Sobel.

El núcleo de suavizado es triangular, con muestras [1 2 1]. Esta es una función triangular, que, muestra, conduce a esos tres valores:

      2 + x ,   if -2 < x < 0
h = { 2     ,   if x = 0
      2 - x ,   if 0 < x < 2

Su derivada es:

            1 ,   if -2 < x < 0
d/dx h = {  0 ,   if x = 0        (not really, but it's the sensible solution)
           -1 ,   if 0 < x < 2

Por lo tanto, podemos ver que la aproximación derivada de la diferencia central puede verse como un muestreo de la derivada analítica de la misma función triangular utilizada para suavizar. Así tenemos:

sob3x3 = [1 2 1]' * d/dx [1 2 1] = d/dx ( [1 2 1]' * [1 2 1] )

Por lo tanto, si desea hacer este núcleo más grande, simplemente agrande el núcleo de suavizado:

sob5x5 = d/dx ( [1 2 3 2 1]' * [1 2 3 2 1] ) = [1 2 3 2 1]' * [1 1 0 -1 -1]

sob7x7 = d/dx ( [1 2 3 4 3 2 1]' * [1 2 3 4 3 2 1] ) = [1 2 3 4 3 2 1]' * [1 1 1 0 -1 -1 -1]

Esto es muy diferente del consejo dado por Adam Bowen , quien sugiere convolviendo el núcleo con el núcleo triangular de 3 pestañas a lo largo de cada dimensión: [1 2 1] * [1 2 1] = [1 4 6 4 1], y [1 2 1] * [1 0 -1] = [1 2 0 -2 -1]. Tenga en cuenta que, debido al teorema del límite central, enrollar este núcleo triangular consigo mismo conduce a un filtro que se aproxima un poco más al gaussiano. Cuanto más grande sea el núcleo que creamos por convoluciones repetidas consigo mismo, más nos aproximamos a este gaussiano. Por lo tanto, en lugar de utilizar este método, también podría muestrear directamente la función gaussiana.

Daniel tiene un largo post en el que sugiere extender el núcleo Sobel de otra manera. La forma del núcleo de suavizado aquí diverge de la aproximación gaussiana, no he tratado de estudiar sus propiedades.

Tenga en cuenta que ninguna de estas tres posibles extensiones del núcleo Sobel son en realidad núcleos Sobel, ya que el núcleo Sobel es explícitamente un núcleo 3x3 ( ver una nota histórica de Sobel sobre su operador, que nunca publicó en realidad).

Obsérvese también que No estoy defendiendo el kernel Sobel extendido derivado aquí. ¡Usa derivados gaussianos!

Gracias por todo, probaré la segunda variante de @Adam Bowen, tome el código C# para Sobel5x5, 7x7, 9x9... matriz generaion para esta variante (tal vez con errores, si encuentra un error o puede optimizar el código - escríbalo allí):

    static void Main(string[] args)
    {
        float[,] Sobel3x3 = new float[,] {
            {-1, 0, 1},
            {-2, 0, 2},
            {-1, 0, 1}};

        float[,] Sobel5x5 = Conv2DforSobelOperator(Sobel3x3);
        float[,] Sobel7x7 = Conv2DforSobelOperator(Sobel5x5);
        Console.ReadKey();
    }

    public static float[,] Conv2DforSobelOperator(float[,] Kernel)
    {
        if (Kernel == null)
            throw new Exception("Kernel = null");

        if (Kernel.GetLength(0) != Kernel.GetLength(1))
            throw new Exception("Kernel matrix must be Square matrix!");

        float[,] BaseMatrix = new float[,] {
            {1, 2, 1},
            {2, 4, 2},
            {1, 2, 1}};

        int KernelSize = Kernel.GetLength(0);
        int HalfKernelSize = KernelSize / 2;
        int OutSize = KernelSize + 2;

        if ((KernelSize & 1) == 0) // Kernel_Size must be: 3, 5, 7, 9 ...
            throw new Exception("Kernel size must be odd (3x3, 5x5, 7x7...)");

        float[,] Out = new float[OutSize, OutSize];
        float[,] InMatrix = new float[OutSize, OutSize];

        for (int x = 0; x < BaseMatrix.GetLength(0); x++)
            for (int y = 0; y < BaseMatrix.GetLength(1); y++)
                InMatrix[HalfKernelSize + x, HalfKernelSize + y] = BaseMatrix[x, y];

        for (int x = 0; x < OutSize; x++)
            for (int y = 0; y < OutSize; y++)
                for (int Kx = 0; Kx < KernelSize; Kx++)
                    for (int Ky = 0; Ky < KernelSize; Ky++)
                    {
                        int X = x + Kx - HalfKernelSize;
                        int Y = y + Ky - HalfKernelSize;

                        if (X >= 0 && Y >= 0 && X < OutSize && Y < OutSize)
                            Out[x, y] += InMatrix[X, Y] * Kernel[KernelSize - 1 - Kx, KernelSize - 1 - Ky];
                    }
        return Out;
    }

Resultados (NormalMap) o se copia allí, donde este metod - №2, @Paul R metod - №1. Ahora estoy usando last, porque da un resultado más suave y es fácil generar núcleos con este código.

Aquí hay una solución simple hecha con python 3 usando numpy y la respuesta @Daniel.

def custom_sobel(shape, axis):
    """
    shape must be odd: eg. (5,5)
    axis is the direction, with 0 to positive x and 1 to positive y
    """
    k = np.zeros(shape)
    p = [(j,i) for j in range(shape[0]) 
           for i in range(shape[1]) 
           if not (i == (shape[1] -1)/2. and j == (shape[0] -1)/2.)]

    for j, i in p:
        j_ = int(j - (shape[0] -1)/2.)
        i_ = int(i - (shape[1] -1)/2.)
        k[j,i] = (i_ if axis==0 else j_)/float(i_*i_ + j_*j_)
    return k

Devuelve el núcleo (5,5) así:

Sobel x:
   [[-0.25 -0.2   0.    0.2   0.25]
    [-0.4  -0.5   0.    0.5   0.4 ]
    [-0.5  -1.    0.    1.    0.5 ]
    [-0.4  -0.5   0.    0.5   0.4 ]
    [-0.25 -0.2   0.    0.2   0.25]]


Sobel y:
   [[-0.25 -0.4  -0.5  -0.4  -0.25]
    [-0.2  -0.5  -1.   -0.5  -0.2 ]
    [ 0.    0.    0.    0.    0.  ]
    [ 0.2   0.5   1.    0.5   0.2 ]
    [ 0.25  0.4   0.5   0.4   0.25]]

Si alguien conoce una mejor manera de hacer eso en python, por favor hágamelo saber. Soy un novato todavía;)