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: [email protected] > To: A Navaei <[email protected]> > > Here is the same example using a compiled function. > > Kent > > > Do you have python-itk-numpy installed? > > > > 2009/2/6 Kent Andre <[email protected]>: > >> > >> 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 > [email protected] > http://www.fenics.org/mailman/listinfo/dolfin-dev
signature.asc
Description: Digital signature
_______________________________________________ DOLFIN-dev mailing list [email protected] http://www.fenics.org/mailman/listinfo/dolfin-dev
