Mentions légales du service

Skip to content
Snippets Groups Projects
Commit 0e13be2d authored by Martin Genet's avatar Martin Genet
Browse files

Cleaning k-space resampling

parent 4d375859
No related branches found
No related tags found
No related merge requests found
......@@ -43,25 +43,28 @@ def downsample_images(
images_ndim = myvtk.getImageDimensionality(
image=image,
verbose=0)
mypy.my_print(verbose, "images_ndim = "+str(images_ndim))
images_spacing = image.GetSpacing()
images_origin = image.GetOrigin()
mypy.my_print(verbose, "images_spacing = "+str(images_spacing))
images_delta = images_spacing[0:images_ndim]
mypy.my_print(verbose, "images_delta = "+str(images_delta))
images_nvoxels = myvtk.getImageDimensions(
image=image,
verbose=0)
mypy.my_print(verbose, "images_nvoxels = "+str(images_nvoxels))
assert (len(images_nvoxels) == images_ndim)
images_npoints = numpy.prod(images_nvoxels)
# mypy.my_print(verbose, "images_npoints = "+str(images_npoints))
assert (images_npoints == image.GetNumberOfPoints())
mypy.my_print(verbose, "images_npoints = "+str(images_npoints))
images_downsampled_nvoxels = numpy.array(images_nvoxels)/numpy.array(downsampling_factors)
images_downsampled_nvoxels = numpy.ceil(images_downsampled_nvoxels)
images_downsampled_nvoxels = [int(n) for n in images_downsampled_nvoxels]
mypy.my_print(verbose, "images_downsampled_nvoxels = "+str(images_downsampled_nvoxels))
downsampling_factors = numpy.array(images_nvoxels)/numpy.array(images_downsampled_nvoxels)
downsampling_factors = list(numpy.array(images_nvoxels)/numpy.array(images_downsampled_nvoxels))
mypy.my_print(verbose, "downsampling_factors = "+str(downsampling_factors))
downsampling_factor = numpy.prod(downsampling_factors)
mypy.my_print(verbose, "downsampling_factor = "+str(downsampling_factor))
images_downsampled_delta = list(numpy.array(images_delta)*numpy.array(downsampling_factors))
# mypy.my_print(verbose, "images_downsampled_delta = "+str(images_downsampled_delta))
images_downsampled_npoints = numpy.prod(images_downsampled_nvoxels)
# mypy.my_print(verbose, "images_downsampled_npoints = "+str(images_downsampled_npoints))
......@@ -149,32 +152,18 @@ def downsample_images(
else:
image_downsampled = vtk.vtkImageData()
if (images_ndim == 1):
extent = [0, images_downsampled_nvoxels[0]-1, 0, 0, 0, 0]
elif (images_ndim == 2):
extent = [0, images_downsampled_nvoxels[0]-1, 0, images_downsampled_nvoxels[1]-1, 0, 0]
elif (images_ndim == 3):
extent = [0, images_downsampled_nvoxels[0]-1, 0, images_downsampled_nvoxels[1]-1, 0, images_downsampled_nvoxels[2]-1]
image_downsampled.SetExtent(extent)
mypy.my_print(verbose, "extent = "+str(extent))
if (images_ndim == 1):
spacing = [images_spacing[0]*downsampling_factors[0], 1., 1.]
elif (images_ndim == 2):
spacing = [images_spacing[0]*downsampling_factors[0], images_spacing[1]*downsampling_factors[1], 1.]
elif (images_ndim == 3):
spacing = [images_spacing[0]*downsampling_factors[0], images_spacing[1]*downsampling_factors[1], images_spacing[2]*downsampling_factors[2]]
image_downsampled.SetSpacing(spacing)
mypy.my_print(verbose, "spacing = "+str(spacing))
if (images_ndim == 1):
origin = [spacing[0]/2, 0., 0.]
elif (images_ndim == 2):
origin = [spacing[0]/2, spacing[1]/2, 0.]
elif (images_ndim == 3):
origin = [spacing[0]/2, spacing[1]/2, spacing[2]/2]
image_downsampled.SetOrigin(origin)
mypy.my_print(verbose, "origin = "+str(origin))
dimensions_downsampled = images_downsampled_nvoxels+[1]*(3-images_ndim)
mypy.my_print(verbose, "dimensions_downsampled = "+str(dimensions_downsampled))
image_downsampled.SetDimensions(dimensions_downsampled)
spacing_downsampled = images_downsampled_delta+[1.]*(3-images_ndim)
mypy.my_print(verbose, "spacing_downsampled = "+str(spacing_downsampled))
image_downsampled.SetSpacing(spacing_downsampled)
origin_downsampled = list(numpy.array(images_downsampled_delta)/2)
origin_downsampled = origin_downsampled+[0.]*(3-images_ndim)
mypy.my_print(verbose, "origin_downsampled = "+str(origin_downsampled))
image_downsampled.SetOrigin(origin_downsampled)
image_downsampled_scalars = myvtk.createDoubleArray(
name="ImageScalars",
......
......@@ -28,29 +28,29 @@ def set_I_woGrad(
image,
X,
I,
vtk_image_scalars,
image_upsampled_scalars,
k_point,
G=None,
Finv=None,
vtk_gradient_vectors=None):
image_gradient_upsampled_vectors=None):
image.I0(X, I)
vtk_image_scalars.SetTuple(k_point, I)
image_upsampled_scalars.SetTuple(k_point, I)
def set_I_wGrad(
image,
X,
I,
vtk_image_scalars,
image_upsampled_scalars,
k_point,
G,
Finv,
vtk_gradient_vectors):
image_gradient_upsampled_vectors):
image.I0_wGrad(X, I, G)
vtk_image_scalars.SetTuple(k_point, I)
image_upsampled_scalars.SetTuple(k_point, I)
G = numpy.dot(G, Finv)
vtk_gradient_vectors.SetTuple(k_point, G)
image_gradient_upsampled_vectors.SetTuple(k_point, G)
def generate_images(
images,
......@@ -90,61 +90,47 @@ def generate_images(
evolution,
generate_image_gradient)
vtk_image = vtk.vtkImageData()
delta = numpy.array(images["L"])/numpy.array(images["n_voxels"])
images["n_voxels_upsampled"] = numpy.array(images["n_voxels"])*numpy.array(images["upsampling_factors"])
delta_upsampled = numpy.array(images["L"])/numpy.array(images["n_voxels_upsampled"])
if (images["n_dim"] == 1):
extent = [0, images["n_voxels_upsampled"][0]-1, 0, 0, 0, 0]
elif (images["n_dim"] == 2):
extent = [0, images["n_voxels_upsampled"][0]-1, 0, images["n_voxels_upsampled"][1]-1, 0, 0]
elif (images["n_dim"] == 3):
extent = [0, images["n_voxels_upsampled"][0]-1, 0, images["n_voxels_upsampled"][1]-1, 0, images["n_voxels_upsampled"][2]-1]
vtk_image.SetExtent(extent)
mypy.my_print(verbose, "extent = "+str(extent))
if (images["n_dim"] == 1):
spacing = [delta_upsampled[0], 1., 1.]
elif (images["n_dim"] == 2):
spacing = [delta_upsampled[0], delta_upsampled[1], 1.]
elif (images["n_dim"] == 3):
spacing = [delta_upsampled[0], delta_upsampled[1], delta_upsampled[2]]
vtk_image.SetSpacing(spacing)
mypy.my_print(verbose, "spacing = "+str(spacing))
if (images["n_dim"] == 1):
origin = [delta[0]/2, 0., 0.]
elif (images["n_dim"] == 2):
origin = [delta[0]/2, delta[1]/2, 0.]
elif (images["n_dim"] == 3):
origin = [delta[0]/2, delta[1]/2, delta[2]/2]
vtk_image.SetOrigin(origin)
mypy.my_print(verbose, "origin = "+str(origin))
n_points = vtk_image.GetNumberOfPoints()
vtk_image_scalars = myvtk.createFloatArray(
image_upsampled = vtk.vtkImageData()
n_voxels_upsampled = list(numpy.array(images["n_voxels"])*numpy.array(images["upsampling_factors"]))
dimensions_upsampled = n_voxels_upsampled+[1]*(3-images["n_dim"])
mypy.my_print(verbose, "dimensions_upsampled = "+str(dimensions_upsampled))
image_upsampled.SetDimensions(dimensions_upsampled)
delta = list(numpy.array(images["L"])/numpy.array(images["n_voxels"]))
delta_upsampled = list(numpy.array(delta)/numpy.array(images["upsampling_factors"]))
spacing_upsampled = delta_upsampled+[1.]*(3-images["n_dim"])
mypy.my_print(verbose, "spacing_upsampled = "+str(spacing_upsampled))
image_upsampled.SetSpacing(spacing_upsampled)
origin_upsampled = list(numpy.array(delta_upsampled)/2)
origin_upsampled = origin_upsampled+[0.]*(3-images["n_dim"])
mypy.my_print(verbose, "origin_upsampled = "+str(origin_upsampled))
image_upsampled.SetOrigin(origin_upsampled)
n_points_upsampled = image_upsampled.GetNumberOfPoints()
image_upsampled_scalars = myvtk.createFloatArray(
name="ImageScalars",
n_components=1,
n_tuples=n_points,
n_tuples=n_points_upsampled,
verbose=verbose-1)
vtk_image.GetPointData().SetScalars(vtk_image_scalars)
image_upsampled.GetPointData().SetScalars(image_upsampled_scalars)
if (generate_image_gradient):
vtk_gradient = vtk.vtkImageData()
vtk_gradient.DeepCopy(vtk_image)
image_gradient_upsampled = vtk.vtkImageData()
image_gradient_upsampled.DeepCopy(image_upsampled)
vtk_gradient_vectors = myvtk.createFloatArray(
image_gradient_upsampled_vectors = myvtk.createFloatArray(
name="ImageScalarsGradient",
n_components=3,
n_tuples=n_points,
n_tuples=n_points_upsampled,
verbose=verbose-1)
vtk_gradient.GetPointData().SetScalars(vtk_gradient_vectors)
image_gradient_upsampled.GetPointData().SetScalars(image_gradient_upsampled_vectors)
else:
vtk_gradient = None
vtk_gradient_vectors = None
image_gradient_upsampled = None
image_gradient_upsampled_vectors = None
x = numpy.empty(3)
X = numpy.empty(3)
......@@ -164,23 +150,23 @@ def generate_images(
t = images["T"]*float(k_frame)/(images["n_frames"]-1) if (images["n_frames"]>1) else 0.
mypy.my_print(verbose, "t = "+str(t))
mapping.init_t(t)
for k_point in xrange(n_points):
vtk_image.GetPoint(k_point, x)
for k_point in xrange(n_points_upsampled):
image_upsampled.GetPoint(k_point, x)
#print "x0 = "+str(x)
mapping.X(x, X, Finv)
#print "X = "+str(X)
set_I(image, X, I, vtk_image_scalars, k_point, G, Finv, vtk_gradient_vectors)
set_I(image, X, I, image_upsampled_scalars, k_point, G, Finv, image_gradient_upsampled_vectors)
global_min = min(global_min, I[0])
global_max = max(global_max, I[0])
#print vtk_image
#print image_upsampled
myvtk.writeImage(
image=vtk_image,
image=image_upsampled,
filename=images["folder"]+"/"+images["basename"]+"_"+str(k_frame).zfill(images["zfill"])+"."+images["ext"],
verbose=verbose-1)
if (generate_image_gradient):
#print vtk_gradient
#print image_gradient_upsampled
myvtk.writeImage(
image=vtk_gradient,
image=image_gradient_upsampled,
filename=images["folder"]+"/"+images["basename"]+"-grad"+"_"+str(k_frame).zfill(images["zfill"])+"."+images["ext"],
verbose=verbose-1)
# mypy.my_print(verbose, "global_min = "+str(global_min))
......@@ -190,6 +176,7 @@ def generate_images(
downsampling = False
else:
downsampling = True
ddic.downsample_images(
images_folder=images["folder"],
images_basename=images["basename"],
......@@ -199,8 +186,10 @@ def generate_images(
verbose=verbose)
if (images["data_type"] in ("float")):
pass
normalizing = False
elif (images["data_type"] in ("unsigned char", "unsigned short", "unsigned int", "unsigned long", "unsigned float", "uint8", "uint16", "uint32", "uint64", "ufloat")):
normalizing = True
ddic.normalize_images(
images_folder=images["folder"],
images_basename=images["basename"]+("_downsampled"*(downsampling)),
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment