Functions | |
| def | new_average_bldr |
| use our new averaging operator | |
| def | OrientationEstimation |
| Takes the image specified by source, and writes the gradient direction angle for the source image (lying between -pi and pi) to the RSF file specified by target. | |
Variables | |
| string | usage = "\nOrientationEstimation.py target=theta.rsf datafile=data.rsf" |
| dictionary | defaults |
| list | __all__ = ['OrientationEstimation'] |
Below, pieces of the actual code are shown in sequence, with a brief description of what the code does preceding the corresponding code snippet.
First, we give the proper usage for this module.
usage="\nOrientationEstimation.py target=theta.rsf datafile=data.rsf"
from SLIMpy import * import numpy as np import rsf as sf from slimproj import *
Set the default value for the input filename containing the initial gradient estimate of the image
defaults = {
'inputgrad_filename':'marmousi_grad.rsf', #Input data as y
}
The new_average_bldr function is a builder function which takes the image to which we wish to apply the sliding average window (z) as its first argument and the size of the averaging window (filtersize) as its second argument.
def new_average_bldr( vectors ):
\"""
use our new averaging operator
\"""
z = vectors[0]
filtersize = vectors[1]
averagefilterx = np.ones((1,filtersize))/filtersize
averagefiltery = np.ones((1,filtersize))/filtersize
Allocate the output (averaged) image
(M,N) = z.shape
x = np.zeros((M,N),'f')
First average along each row, and then along each column. This is equivalent to taking the average in a rectangular (in this case square) window across the image.
for i in range(0,M):
x[i,:] = np.convolve(z[i,:],np.ones(filtersize)/filtersize,mode='same')
for j in range(0,N):
x[:,j] = np.convolve(x[:,j],np.ones(filtersize)/filtersize,mode='same')
Return the image computed by taking windowed averages.
return x
Read in input image.
def OrientationEstimation( target, source, env ):
inputfilename = source[0]
inputimage = sf.Input(source[0])
M is the number of rows and N is the number of columns in image
M = inputimage.int("n2")
N = inputimage.int("n1")
Read in previously computed partial derivative in x direction
Jx = np.zeros((M,N),'f')
noisyinputgradx_filename = source[1]
noisyinputgradxrsf = sf.Input(noisyinputgradx_filename)
noisyinputgradxrsf.read(Jx)
Read in previously computed partial derivative in y direction
Jy = np.zeros((M,N),'f')
noisyinputgrady_filename = source[2]
noisyinputgradyrsf = sf.Input(noisyinputgrady_filename)
noisyinputgradyrsf.read(Jy)
Allocate space for elements of structure tensor and initialize to zero.
Jtensorxx = np.zeros((M,N))
Jtensorxy = np.zeros((M,N))
Jtensoryy = np.zeros((M,N))
Form elements of structure tensor (Outer product of gradient vector with itself).
Jtensorxx = Jx*Jx
Jtensorxy = Jx*Jy
Jtensoryy = Jy*Jy
Replace each pixel's structure tensor with a local average of the structure tensor around the pixel.
averagefilterwidth = 5
Jtensorxx = new_average_bldr([Jtensorxx,averagefilterwidth])
Jtensorxy = new_average_bldr([Jtensorxy,averagefilterwidth])
Jtensoryy = new_average_bldr([Jtensoryy,averagefilterwidth])
Estimated gradient orientation is in the same direction as the eigenvector corresponding to the largest eigenvalue of the structure tensor.
Jtensoreig = np.zeros((2,M,N),'f')
Jmag = np.zeros((M,N),'f')
tempJ = np.zeros((2,2),'f')
for i in range(0,M):
for j in range(0,N):
tempJ[0,0] = Jtensorxx[i,j]
tempJ[0,1] = Jtensorxy[i,j]
tempJ[1,0] = Jtensorxy[i,j]
tempJ[1,1] = Jtensoryy[i,j]
tempJeigs = np.linalg.eig(tempJ)
if (tempJeigs[0][0] > tempJeigs[0][1]):
Jtensoreig[0][i][j] = tempJeigs[1][0][0]
Jtensoreig[1][i][j] = tempJeigs[1][0][1]
else:
Jtensoreig[0][i][j] = tempJeigs[1][1][0]
Jtensoreig[1][i][j] = tempJeigs[1][1][1]
Jmag[i][j] = np.sqrt(tempJeigs[0][0]+tempJeigs[0][1])
Compute the gradient angle from the elements of the gradient vector.
theta = np.array((M,N),'f')
theta = np.arctan(Jtensoreig[1,:,:]/Jtensoreig[0,:,:])+np.pi/2
Write the estimated gradient orientation to the RSF File specified by target.
output = sf.Output(target)
output.put("n1",N)
output.put("n2",M)
output.write(theta)
End the current function.
End()
Finally, export the OrientationEstimation function so that it can be used by other programs.
__all__ = ['OrientationEstimation']
from slimproj_core.builders.CreateBuilders import add_to_slim_env
add_to_slim_env("OrientationEstimation",OrientationEstimation)
| def OrientationEstimation.new_average_bldr | ( | vectors | ) |
| def OrientationEstimation.OrientationEstimation | ( | target, | ||
| source, | ||||
| env | ||||
| ) |
Takes the image specified by source, and writes the gradient direction angle for the source image (lying between -pi and pi) to the RSF file specified by target.
Definition at line 241 of file OrientationEstimation.py.
| list __all__ = ['OrientationEstimation'] |
Definition at line 323 of file OrientationEstimation.py.
| dictionary defaults |
Initial value:
{
'inputgrad_filename':'marmousi_grad.rsf', #Input data as y
}
Definition at line 212 of file OrientationEstimation.py.
| string usage = "\nOrientationEstimation.py target=theta.rsf datafile=data.rsf" |
Definition at line 204 of file OrientationEstimation.py.
1.5.6