dolfin team mailing list archive

[Anders Logg <logg@xxxxxxxxx>]

How about adding ImageFunction to SpecialFunctions.h/cpp.

It would require using the ITK C++ interface (I assume there is one),
adding a check for ITK to the build system and a flag HAS_ITK.

-- 
Anders


On Mon, Feb 09, 2009 at 12:33:14PM +0100, Kent Andre wrote:
> 

Content-Description: Forwarded message - Re: [DOLFIN-dev] Image to Function data structure conversion
> Date: Sat, 7 Feb 2009 20:52:38 +0100 (CET)
> Subject: Re: [DOLFIN-dev] Image to Function data structure conversion
> From: kent-and@xxxxxxxxxxxx
> To: A Navaei <axnavaei@xxxxxxxxxxxxxx>
> 
> 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()
> 
> 
> 


> _______________________________________________
> DOLFIN-dev mailing list
> DOLFIN-dev@xxxxxxxxxx
> http://www.fenics.org/mailman/listinfo/dolfin-dev

Attachment: signature.asc
Description: Digital signature