Comment une convolution peut-elle être exprimée comme une multiplication matricielle (forme matricielle)?

11

Je sais que cette question peut ne pas être très pertinente pour la programmation, mais si je ne comprends pas la théorie derrière le traitement d'image, je ne pourrai jamais mettre en œuvre quelque chose dans la pratique.

Si je comprends bien, les filtres gaussiens sont convolués avec une image pour la réduction du bruit car ils calculent une moyenne pondérée du voisinage d'un pixel et ils sont très utiles pour la détection des contours, car vous pouvez appliquer un flou et dériver l'image en même temps en convoluant simplement avec la dérivée d'une fonction gaussienne.

Mais quelqu'un peut-il m'expliquer ou me donner des références sur la façon dont ils sont calculés?

Par exemple, le détecteur de bord de Canny parle d'un filtre gaussien 5x5, mais comment ont-ils obtenu ces nombres particuliers? Et comment sont-ils passés d'une convolution continue à une multiplication matricielle?

Matteo
la source
J'ai ajouté une réponse avec le code complet pour générer une matrice pour la convolution d'image.
Royi

Réponses:

3

Pour que cette opération fonctionne, vous devez imaginer que votre image est remodelée en tant que vecteur. Ensuite, ce vecteur est multiplié à sa gauche par la matrice de convolution afin d'obtenir l'image floue. Notez que le résultat est également un vecteur de la même taille que l'entrée, c'est-à-dire une image de la même taille.

Chaque ligne de la matrice de convolution correspond à un pixel dans l'image d'entrée. Il contient le poids des contributions de tous les autres pixels de l'image à la contrepartie floue du pixel considéré.

Prenons un exemple: boîte flou de taille pixels sur une image de taille 6 × 6 pixels. L'image remodelée est une colonne de 36 choix, tandis que la matrice de flou a une taille 36 × 36 .3×36×636×36

  • Initions cette matrice à 0 partout.
  • Maintenant, considérez le pixel de coordonnées dans l'image d'entrée (pas sur sa bordure pour plus de simplicité). Son homologue floue est obtenue en appliquant un poids de 1 / neuf à lui - même et chacun de ses voisins au niveau des positions ( i - 1 , j - 1 ) ; ( i - 1 , j ) , ( i - 1 , j + 1 ) , , ( i + 1 ,(je,j)1/9 .(je-1,j-1);(je-1,j),(je-1,j+1),,(je+1,j+1)
  • Dans le vecteur colonne, le pixel a la position 6 i + j (en supposant un ordre lexicographique). Nous rapportons le poids une / 9 dans le ( 6 i + j ) -ième ligne de la matrice de flou.(je,j)6je+j1/9(6je+j)
  • Faites de même avec tous les autres pixels.

Une illustration visuelle d'un processus étroitement lié (convolution + soustraction) peut être trouvée sur ce billet de blog (à partir de mon blog personnel).

sansuiso
la source
un lien est mort.
gauteh
2

Pour les applications aux images ou aux réseaux de convolution, pour utiliser plus efficacement les multiplicateurs matriciels dans les GPU modernes, les entrées sont généralement remodelées en colonnes d'une matrice d'activation qui peut ensuite être multipliée avec plusieurs filtres / noyaux à la fois.

Consultez ce lien depuis le CS231n de Stanford et faites défiler jusqu'à la section «Implémentation en tant que multiplication matricielle» pour plus de détails.

Le processus fonctionne en prenant tous les correctifs locaux sur une image d'entrée ou une carte d'activation, ceux qui seraient multipliés avec le noyau, et en les étirant dans une colonne d'une nouvelle matrice X via une opération communément appelée im2col. Les noyaux sont également étirés pour remplir les rangées d'une matrice de poids W de sorte que lors de l'exécution de l'opération de matrice W * X, la matrice résultante Y ait tous les résultats de la convolution. Enfin, la matrice Y doit être remodelée à nouveau en reconvertissant les colonnes en images par une opération généralement appelée cal2im.

Rodrigo
la source
1
Ceci est un très bon lien, merci! Cependant, il est recommandé d'ajouter les extraits importants du lien dans la réponse, de cette façon la réponse est valide même si le lien est rompu. Veuillez envisager de modifier votre réponse pour la faire accepter!
Matteo
1

La convolution dans le domaine temporel est égale à la multiplication matricielle dans le domaine fréquentiel et vice versa.

Le filtrage est équivalent à la convolution dans le domaine temporel et donc à la multiplication matricielle dans le domaine fréquentiel.

Quant aux cartes ou masques 5x5, ils proviennent de la discrétisation des opérateurs canny / sobel.

Naresh
la source
2
Je ne suis pas d'accord avec le fait que le filtrage est une convolution dans le domaine fréquentiel. Le type de filtres dont nous parlons ici sont des convolutions dans le domaine spatial (c'est-à-dire la multiplication par élément par la réponse du filtre dans le domaine fréquentiel).
pichenettes
Merci d'avoir corrigé ce que j'ai écrit. J'ai fait un montage ultérieur. Je suppose que je devrais vérifier mes réponses avant de poster. Cependant, la majorité de ma réponse est toujours valable.
Naresh
La transformée de Fourier transforme en effet les convolutions en multiplications (et vice versa). Cependant, ce sont des multiplications par pinte, alors que la question concerne les multiplications matricielles-vectorielles qui sont obtenues en remodelant les images.
sansuiso
J'ai mentionné à quel point la discrétisation des opérateurs est à l'origine des matrices 5x5 obtenues pour les opérateurs canny / sobel.
Naresh
1

J'ai écrit une fonction qui résout cela dans mon référentiel GitHub StackOverflow Q2080835 (Jetez un oeil à CreateImageConvMtx()).
En fait, la fonction peut prendre en charge n'importe quelle forme de convolution que vous souhaitez - full, sameet valid.

Le code est le suivant:

function [ mK ] = CreateImageConvMtx( mH, numRows, numCols, convShape )

CONVOLUTION_SHAPE_FULL  = 1;
CONVOLUTION_SHAPE_SAME  = 2;
CONVOLUTION_SHAPE_VALID = 3;

switch(convShape)
    case(CONVOLUTION_SHAPE_FULL)
        % Code for the 'full' case
        convShapeString = 'full';
    case(CONVOLUTION_SHAPE_SAME)
        % Code for the 'same' case
        convShapeString = 'same';
    case(CONVOLUTION_SHAPE_VALID)
        % Code for the 'valid' case
        convShapeString = 'valid';
end

mImpulse = zeros(numRows, numCols);

for ii = numel(mImpulse):-1:1
    mImpulse(ii)    = 1; %<! Create impulse image corresponding to i-th output matrix column
    mTmp            = sparse(conv2(mImpulse, mH, convShapeString)); %<! The impulse response
    cColumn{ii}     = mTmp(:);
    mImpulse(ii)    = 0;
end

mK = cell2mat(cColumn);


end

J'ai également créé une fonction pour créer une matrice pour le filtrage d'images (idées similaires à MATLAB imfilter()):

function [ mK ] = CreateImageFilterMtx( mH, numRows, numCols, operationMode, boundaryMode )
%UNTITLED6 Summary of this function goes here
%   Detailed explanation goes here

OPERATION_MODE_CONVOLUTION = 1;
OPERATION_MODE_CORRELATION = 2;

BOUNDARY_MODE_ZEROS         = 1;
BOUNDARY_MODE_SYMMETRIC     = 2;
BOUNDARY_MODE_REPLICATE     = 3;
BOUNDARY_MODE_CIRCULAR      = 4;

switch(operationMode)
    case(OPERATION_MODE_CONVOLUTION)
        mH = mH(end:-1:1, end:-1:1);
    case(OPERATION_MODE_CORRELATION)
        % mH = mH; %<! Default Code is correlation
end

switch(boundaryMode)
    case(BOUNDARY_MODE_ZEROS)
        mK = CreateConvMtxZeros(mH, numRows, numCols);
    case(BOUNDARY_MODE_SYMMETRIC)
        mK = CreateConvMtxSymmetric(mH, numRows, numCols);
    case(BOUNDARY_MODE_REPLICATE)
        mK = CreateConvMtxReplicate(mH, numRows, numCols);
    case(BOUNDARY_MODE_CIRCULAR)
        mK = CreateConvMtxCircular(mH, numRows, numCols);
end


end


function [ mK ] = CreateConvMtxZeros( mH, numRows, numCols )
%UNTITLED6 Summary of this function goes here
%   Detailed explanation goes here

numElementsImage    = numRows * numCols;
numRowsKernel       = size(mH, 1);
numColsKernel       = size(mH, 2);
numElementsKernel   = numRowsKernel * numColsKernel;

vRows = reshape(repmat(1:numElementsImage, numElementsKernel, 1), numElementsImage * numElementsKernel, 1);
vCols = zeros(numElementsImage * numElementsKernel, 1);
vVals = zeros(numElementsImage * numElementsKernel, 1);

kernelRadiusV = floor(numRowsKernel / 2);
kernelRadiusH = floor(numColsKernel / 2);

pxIdx       = 0;
elmntIdx    = 0;

for jj = 1:numCols
    for ii = 1:numRows
        pxIdx = pxIdx + 1;
        for ll = -kernelRadiusH:kernelRadiusH
            for kk = -kernelRadiusV:kernelRadiusV
                elmntIdx = elmntIdx + 1;

                pxShift = (ll * numCols) + kk;

                if((ii + kk <= numRows) && (ii + kk >= 1) && (jj + ll <= numCols) && (jj + ll >= 1))
                    vCols(elmntIdx) = pxIdx + pxShift;
                    vVals(elmntIdx) = mH(kk + kernelRadiusV + 1, ll + kernelRadiusH + 1);
                else
                    vCols(elmntIdx) = pxIdx;
                    vVals(elmntIdx) = 0; % See the accumulation property of 'sparse()'.
                end
            end
        end
    end
end

mK = sparse(vRows, vCols, vVals, numElementsImage, numElementsImage);


end


function [ mK ] = CreateConvMtxSymmetric( mH, numRows, numCols )
%UNTITLED6 Summary of this function goes here
%   Detailed explanation goes here

numElementsImage    = numRows * numCols;
numRowsKernel       = size(mH, 1);
numColsKernel       = size(mH, 2);
numElementsKernel   = numRowsKernel * numColsKernel;

vRows = reshape(repmat(1:numElementsImage, numElementsKernel, 1), numElementsImage * numElementsKernel, 1);
vCols = zeros(numElementsImage * numElementsKernel, 1);
vVals = zeros(numElementsImage * numElementsKernel, 1);

kernelRadiusV = floor(numRowsKernel / 2);
kernelRadiusH = floor(numColsKernel / 2);

pxIdx       = 0;
elmntIdx    = 0;

for jj = 1:numCols
    for ii = 1:numRows
        pxIdx = pxIdx + 1;
        for ll = -kernelRadiusH:kernelRadiusH
            for kk = -kernelRadiusV:kernelRadiusV
                elmntIdx = elmntIdx + 1;

                pxShift = (ll * numCols) + kk;

                if(ii + kk > numRows)
                    pxShift = pxShift - (2 * (ii + kk - numRows) - 1);
                end

                if(ii + kk < 1)
                    pxShift = pxShift + (2 * (1 -(ii + kk)) - 1);
                end

                if(jj + ll > numCols)
                    pxShift = pxShift - ((2 * (jj + ll - numCols) - 1) * numCols);
                end

                if(jj + ll < 1)
                    pxShift = pxShift + ((2 * (1 - (jj + ll)) - 1) * numCols);
                end

                vCols(elmntIdx) = pxIdx + pxShift;
                vVals(elmntIdx) = mH(kk + kernelRadiusV + 1, ll + kernelRadiusH + 1);

            end
        end
    end
end

mK = sparse(vRows, vCols, vVals, numElementsImage, numElementsImage);


end


function [ mK ] = CreateConvMtxReplicate( mH, numRows, numCols )
%UNTITLED6 Summary of this function goes here
%   Detailed explanation goes here

numElementsImage    = numRows * numCols;
numRowsKernel       = size(mH, 1);
numColsKernel       = size(mH, 2);
numElementsKernel   = numRowsKernel * numColsKernel;

vRows = reshape(repmat(1:numElementsImage, numElementsKernel, 1), numElementsImage * numElementsKernel, 1);
vCols = zeros(numElementsImage * numElementsKernel, 1);
vVals = zeros(numElementsImage * numElementsKernel, 1);

kernelRadiusV = floor(numRowsKernel / 2);
kernelRadiusH = floor(numColsKernel / 2);

pxIdx       = 0;
elmntIdx    = 0;

for jj = 1:numCols
    for ii = 1:numRows
        pxIdx = pxIdx + 1;
        for ll = -kernelRadiusH:kernelRadiusH
            for kk = -kernelRadiusV:kernelRadiusV
                elmntIdx = elmntIdx + 1;

                pxShift = (ll * numCols) + kk;

                if(ii + kk > numRows)
                    pxShift = pxShift - (ii + kk - numRows);
                end

                if(ii + kk < 1)
                    pxShift = pxShift + (1 -(ii + kk));
                end

                if(jj + ll > numCols)
                    pxShift = pxShift - ((jj + ll - numCols) * numCols);
                end

                if(jj + ll < 1)
                    pxShift = pxShift + ((1 - (jj + ll)) * numCols);
                end

                vCols(elmntIdx) = pxIdx + pxShift;
                vVals(elmntIdx) = mH(kk + kernelRadiusV + 1, ll + kernelRadiusH + 1);

            end
        end
    end
end

mK = sparse(vRows, vCols, vVals, numElementsImage, numElementsImage);


end


function [ mK ] = CreateConvMtxCircular( mH, numRows, numCols )
%UNTITLED6 Summary of this function goes here
%   Detailed explanation goes here

numElementsImage    = numRows * numCols;
numRowsKernel       = size(mH, 1);
numColsKernel       = size(mH, 2);
numElementsKernel   = numRowsKernel * numColsKernel;

vRows = reshape(repmat(1:numElementsImage, numElementsKernel, 1), numElementsImage * numElementsKernel, 1);
vCols = zeros(numElementsImage * numElementsKernel, 1);
vVals = zeros(numElementsImage * numElementsKernel, 1);

kernelRadiusV = floor(numRowsKernel / 2);
kernelRadiusH = floor(numColsKernel / 2);

pxIdx       = 0;
elmntIdx    = 0;

for jj = 1:numCols
    for ii = 1:numRows
        pxIdx = pxIdx + 1;
        for ll = -kernelRadiusH:kernelRadiusH
            for kk = -kernelRadiusV:kernelRadiusV
                elmntIdx = elmntIdx + 1;

                pxShift = (ll * numCols) + kk;

                if(ii + kk > numRows)
                    pxShift = pxShift - numRows;
                end

                if(ii + kk < 1)
                    pxShift = pxShift + numRows;
                end

                if(jj + ll > numCols)
                    pxShift = pxShift - (numCols * numCols);
                end

                if(jj + ll < 1)
                    pxShift = pxShift + (numCols * numCols);
                end

                vCols(elmntIdx) = pxIdx + pxShift;
                vVals(elmntIdx) = mH(kk + kernelRadiusV + 1, ll + kernelRadiusH + 1);

            end
        end
    end
end

mK = sparse(vRows, vCols, vVals, numElementsImage, numElementsImage);


end

Le code a été validé par rapport à MATLAB imfilter().

Le code complet est disponible dans mon référentiel GitHub StackOverflow Q2080835 .

Royi
la source