[vtkusers] segfault with MarchingCubes

[David Gobbi]

Hi Clément,

The "number of scalar components" should be 1 for your data:

  dataImporter.SetNumberOfScalarComponents(1)

There are a couple places where you can improve efficiency:

  # directly create an array of the desired type
  npa=np.zeros((150,150,150), np.uint8)

  # directly pass the array to VTK (without tostring())
  dataImporter.CopyImportVoidPointer(npa, npa.size)

Cheers,
 - David


On Mon, Jul 10, 2017 at 2:29 AM, Clément D. <clement.douarre at gmail.com>
wrote:

> Hello all,
>
> I am trying to use the Marching Cubes algorithm to show a 3D surface of a
> 3D image given by a Python (numpy array).
> I thought I had it all mapped out, things were going great with arrays of
> small size (generated and rendered 3D isosurfaces), but now that I'm
> working with arrays of slightly bigger size, the programm crashes with
> Python error 11 SIGSEGV which means seg fault. What worries me a bit is
> that the array I'm using is not that big, 150*150*150, and I have a fairly
> recent computer who shouldn't have any problem handling it.
>
> Code MWE is below (in Python)
>
> =========================================== Code snippet
> =======================
>
> import numpy as np
> import vtk
>
> # Creating a numpy array representing a 3D image
> npa=np.zeros((150,150,150)).astype(np.uint8)
> npa[100:102,100:102,100:102]=255
>
> # Importing it as  object vtk (based on vtk example
> http://www.vtk.org/Wiki/VTK/Examples/Python/vtkWithNumpy)
> dataImporter = vtk.vtkImageImport()
> dataString = npa.tostring()
> [hh, ww, dd] = npa.shape
> dataImporter.CopyImportVoidPointer(dataString, len(dataString))
> dataImporter.SetDataScalarTypeToUnsignedChar()
> dataImporter.SetNumberOfScalarComponents(3)
> dataImporter.SetDataExtent(0, dd - 1, 0, ww - 1, 0, hh - 1)
> dataImporter.SetWholeExtent(0, dd - 1, 0, ww - 1, 0, hh - 1)
> dataImporter.Update()
> vol = dataImporter.GetOutputPort()
> # (print(dataImporter.GetOutput().GetClassName()) yields "vtk.ImageData")
>
> # Use MarchingCubes to generate iso-surface
> surface = vtk.vtkMarchingCubes();
> surface.SetInputConnection(vol);
> surface.ComputeNormalsOn();
> surface.ComputeScalarsOn();
> isoValue = 128
> surface.SetValue(0, isoValue);
>
> #This line causes segfault if npa is "too big"
> surface.Update();
>
> # .... rest of pipeline
> # decimator = vtk.vtkDecimatePro(); ....
> # smoother = vtk.vtkSmoothPolyDataFilter(); ....
> # connectivityFilter = vtk.vtkPolyDataConnectivityFilter(); ....
> # mesh = vtk.vtkPolyData(); ....
> # mapper = vtk.vtkPolyDataMapper(); ....
> # actor = vtk.vtkActor(); ....
>
>
> =========================================== End of code snippet
> =======================
>
>
>
> So I have a couple of questions here :
> 1/ Is this a problem of input for MarchingCubes being too big ? Or is this
> segfault indicating something else ?
> - arguments for the former : Reducing the size of npa in this code makes
> it work
> - arguments for the latter : This code I'm showing is a MWE for my code
> which is more complicated, and in my code, this code portion fails with a
> large array, but succeeds if I reduce in size *another* array which has
> nothing to do with MarchingCubes() (again, that other array is not that
> big, so I don't think I'm reaching the limits of Python memory or anything
> like that...)
>
> 2/ Am I doing the whole import numpy thing correctly ? I've based myself
> on the vtk example, but I am unsure of whether it is the recommended way
> for my case or if it was case specific
>
>
> Thanks everyone,
> Clément
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://public.kitware.com/pipermail/vtkusers/attachments/20170710/cd751ba9/attachment.html>