--- Begin Message ---
Here is the same example using a compiled function.
Kent
> Do you have python-itk-numpy installed?
>
> 2009/2/6 Kent Andre <kent-and@xxxxxxxxx>:
>>
>> Great, but I get the following problem:
>>
>>
>>
>> Traceback (most recent call last):
>> File "itk-dolfin.py", line 16, in <module>
>> itk2numpy = itk.PyBuffer[inType]
>> File "/usr/lib/InsightToolkit/WrapITK/Python/itkLazy.py", line 14, in
>> __getattribute__
>> value = types.ModuleType.__getattribute__(self, attr)
>> AttributeError: 'LazyITKModule' object has no attribute 'PyBuffer'
>>
>> On to., 2009-02-05 at 19:18 +0000, A Navaei wrote:
>>> >
>>> > Let me know when you have created debian packages. It seems I have to
>>> > rebuild everything, even cmake. Compilation is in principle no
>>> problem,
>>> > but I already compile a decent number of packages on a regular basis
>>> and
>>> > being non-bleeding edge on packages I don't know well seems tempting.
>>> >
>>> > I might download and compile, but not today.
>>> >
>>> > Kent
>>> >
>>> >
>>>
>>> python-itk debian, based on the new WrapITK, is just published:
>>> http://code.google.com/p/wrapitk/wiki/WrapITKBinaries .The 32 bit
>>> version is ready, the 64 bit version is being built on opensuse build
>>> service and will be available in a few hours.
>>>
>>> I have tried the [itk]->[numpy]->[dolfin] example, and it already
>>> looks great, here is the output:
>>>
>>> Checking mesh ordering (finished).
>>> 211 * 212 = 44732 pixels took 246.780 ms transfering from itk to
>>> dolfin through numpy
>>> That is 5.516854 us per pixel.
>>>
>>> I am going to add the itk-dolfin interface in c++ in the next few
>>> days, which should make the conversion even faster leading to more
>>> reasonable performance for 3D data. Below is an updated version of the
>>> example. Let me know if there are any problems.
>>>
>>>
>>> -Ali
>>>
>>> ----------------------------------------------------------
>>> # [itk]->[numpy]->[dolfin] test
>>>
>>> from dolfin import *
>>> from numpy import *
>>> import itk
>>> import time
>>>
>>> dim = 2
>>> inType = itk.Image[itk.D, dim]
>>>
>>> reader = itk.ImageFileReader[inType].New(FileName="/path/to/file.bmp")
>>> reader.Update()
>>>
>>> t1 = time.time()
>>>
>>> itk2numpy = itk.PyBuffer[inType]
>>> numpy_arr = itk2numpy.GetArrayFromImage( reader.GetOutput() )
>>> shape = numpy_arr.shape
>>>
>>> mesh = UnitSquare(shape[0]-1, shape[1]-1)
>>>
>>> class ImageFunction(Function):
>>> def eval(self, value, x):
>>> i = int((self.A.shape[0]-1)*x[0])
>>> j = int((self.A.shape[1]-1)*x[1])
>>> # print i,j
>>> value[0] = self.A[(i,j)]
>>>
>>> V = FunctionSpace(mesh, "CG", 1)
>>> f = ImageFunction(V)
>>> f.A = numpy_arr
>>>
>>> t2 = time.time()
>>> print '%i * %i = %i pixels took %0.3f ms transfering from itk to
>>> dolfin through numpy'% (shape[0], shape[1], shape[0]*shape[1],
>>> (t2-t1)*1e3)
>>> print 'That is %f us per pixel.'% ((t2-t1) * 1e6 / (shape[0]*shape[1]))
>>>
>>> plot(f) # can viper plot 2d images better than this?
>>> interactive()
>>
>>
>
from dolfin import *
from numpy import *
import itk
import time
dim = 2
inType = itk.Image[itk.D, dim]
reader = itk.ImageFileReader[inType].New(FileName="bowling.jpg")
reader.Update()
t1 = time.time()
itk2numpy = itk.PyBuffer[inType]
numpy_arr = itk2numpy.GetArrayFromImage( reader.GetOutput() )
shape = numpy_arr.shape
mesh = UnitSquare(shape[0]-1, shape[1]-1)
code = """
class ItkFunc : public Function
{
private:
int N;
int M;
double *image;
public:
void set_image(PyObject* a_){
// Do some more type checks ...
PyArrayObject* a=(PyArrayObject*) a_;
N = a->dimensions[0];
M = a->dimensions[1];
image = (double*) a->data;
}
ItkFunc(FunctionSpace& V) : Function(V), N(0), M(0) {}
void eval(double* values, const double* x) const
{
int i = int((N-1)*x[0]);
int j = int((M-1)*x[1]);
int ij = i + N*j;
values[0] = image[i + (N-1)*j];
}
};
"""
V = FunctionSpace(mesh, "CG", 1)
f = Function(V, code)
f.set_image(numpy_arr)
t2 = time.time()
print '%i * %i = %i pixels took %0.3f ms transfering from itk to dolfin through numpy'% (shape[0], shape[1], shape[0]*shape[1], (t2-t1)*1e3)
print 'That is %f us per pixel.'% ((t2-t1) * 1e6 / (shape[0]*shape[1]))
plot(f) # can viper plot 2d images better than this?
interactive()
--- End Message ---
