I've been messing around with Matlab -- well, Octave really -- trying to reverse engineer exactly how
Dr. Neil Krawetz creates his "luminance gradient" analysis images like the one in
this blog entry (the second of the three analysis images). This is not an easy task, because I don't know Matlab's syntax or features whatsoever, I have only vague and limited linear algebra knowledge from a good 15+ years ago, and I have only some very basic knowledge of digital signal processing. Fairly fundamental tools like Fourier transforms are still beyond my reach at the moment, for instance.
But I thought, and tinkered, and thought some more, and I eventually stumbled into a solution that was very close but not quite right. I computed the luminance of each pixel by simple average: (R + G + B)/3. Then I subsampled the resulting grayscale image at 3x3. I figured I could compute a vector the center pixel by treating each of the surrounding pixels as a vector (i.e., an "arrow") pointing outward, each with a length equal to their luminance. All I would need to do is separate out the X and Y components with sin/cos and sum them into the resulting vector. So if this is a 3x3 subsample:
Code:
D C B
E . A
F G H
The luminance vector should be:
Code:
Lx=cosd(0)*A + cosd(45)*B + cosd(90)*C + cosd(135)*D + cosd(180)*E + cosd(225)*F + cosd(270)*G + cosd(315)*H
Ly=sind(0)*A + sind(45)*B + sind(90)*C + sind(135)*D + sind(180)*E + sind(225)*F + sind(270)*G + sind(315)*H
Which simplifies down to:
Code:
Lx=A + 0.7071B - 0.7071D - E-0.7071F + 0.7071H
Ly=0.7071B + C + 0.7071D - 0.7071F - G - 0.7071H
Doing this processing with nested
for loops in Matlab came up with some not-quite-right-but-close images, but was absolutely slow as balls because of the way that Matlab scripts work (purely interpreted with no JIT compiling, or even loop caching). It would take a couple minutes for it to process just a single 800x600-ish image. So I started poking around in GIMP, trying to see if I could come up with a way to do this operation faster. When I ran into the convolution matrix feature, I realized I had a winner. Lx and Ly were just:
Code:
[ -0.7071 0 0.7071 [ 0.7071 1 0.7071
-1 0 1 and 0 0 0
-0.7071 0 0.7071 ] -0.7071 -1 -0.7071 ]
Googling a bit more, I found that Matlab could also do convolution matrices via the imfilter() function. Ah, hah! No more ridiculously slow for-loops. It was really fast now, but the results still weren't quite the same as Krawetz's. Then I noticed that Matlab's help page for imfilter referenced the fspecial() function, which would prepare various stock filter matrices for you to plug in to the imfilter() command. I looked up the help on fspecial() and sifted through the options. It didn't take long to notice two pre-defined filters called "sobel" and "prewitt".
Sobel:
Code:
[ 1 2 1
0 0 0
-1 -2 -1 ]
Prewitt:
Code:
[ 1 1 1
0 0 0
-1 -1 -1]
Sonofa! Yeah, so I just spent 3 or 4 hours more or less reinventing the
Sobel and
Prewitt operators.
:facepalm:
I don't know whether this should make me feel smart or dumb.
Login
Register
FAQ
Search
