2

If I have a contour in Matlab obtained

 [f, v] = isosurface(x, y, z, v, isovalue)

is there a clean way to apply a transformation to the surfaces and nicely plot the result as a smooth surface? The transformation T is nonlinear.

I tried to apply the transformation T to both f and vert and use patch but couldn't quite get it to work.

Amro
  • 123,847
  • 25
  • 243
  • 454
db1234
  • 797
  • 2
  • 11
  • 28
  • 1
    You should explain more about what exactly you want to do. What is the T transformation? – Adiel Sep 18 '14 at 00:28
  • Assume that T is any general transformation taking points in 3D space to points in 3D space, e.g. T(x,y,z) = (xy, z^2 cos(x), e^y), though my specific is more complicated and not worth describing. I want to take the surface S extracted from isosurface, apply T to S, yielding a new surface T(S). I want a nice way to plot T(S) as a smooth surface. – db1234 Sep 18 '14 at 13:05
  • 3
    I don't think that would work (applying a transformation on the resulting iso-surface vertices/faces). Perhaps you can apply the transformation on the volume data first, and then plot the [isosurface](https://en.wikipedia.org/wiki/Isosurface) of the transformed data? – Amro Sep 20 '14 at 15:06
  • I am not sure what you mean by "apply the transformation on the volume data first". By volume data, I assume you mean "V", which is a scalar field, and I cannot apply a transformation to it. – db1234 Sep 22 '14 at 08:13
  • I think he meant applying the transformation `T` on all the initial volume coordinates `x,y,z` (`[X,Y,Z]=T(x,y,z)`, then pull the isosurface on these new coordinates (`[f, v] = isosurface(X, Y, Z, v, isovalue)`). Regarding the value of the field `v` it may or may not make sense to apply the transformation to it (if even possible), that depends only on your specific problem. – Hoki Sep 22 '14 at 10:33
  • @Hoki: that's what I was thinking of, yes. In that approach we might need to use something like `griddatan` to interpolate the volume from the new transformed X/Y/Z coordinates, since `isosurface` requires uniformly-gridded data as far as I know. Now I'm not sure that will give the same results as what you've shown in your answer, where you simply transformed the isosurface vertices and kept the connecting faces... – Amro Sep 22 '14 at 14:10
  • Recall that an isosurface represents a level set of a continuous 3D function `f(x,y,z)=v` where the function has a specified constant value. So by moving the vertices around with the faces fixed, the new surface no longer represents the expected level set... – Amro Sep 22 '14 at 14:11
  • @Amro, I was curious about that too. I tried it. I applied the transformation on the original data then queried the isosurface (for the same level) and I get exactly the same result. I guess it makes sense mathematically. Now I have no idea if pulling the same level of isosurface on transformed coordinates makes sense in "physics" term (that depends on the problem of the OP) ... But if it does, then it is much faster to transform only the isosurface (in my example it's only 3x3105 points to transform versus the 3x31250 points of the full volume) – Hoki Sep 23 '14 at 11:04

1 Answers1

8

The trick is to apply the transformation on your vertices, but keep the same faces data. This way the faces always link the same points, regardless of their new positions.

Since there are no sample data I took the Matlab example as a starting point. This is coming from the Matlab isosurface page (very slightly modified for this example):

%// Generate an isosurface
[x,y,z,v] = flow;
fv = isosurface(x,y,z,v,-3) ;
figure(1);cla
p1 = patch(fv,'FaceColor','red','EdgeColor','none');
%// refine the view
grid off ; set(gca,'Color','none') ; daspect([1,1,1]) ; view(3) ; axis tight ; camlight ; lighting gouraud

This output:
Original isosurface

Nothing original so far. Just note that I used the single structure output type fv instead of the 2 separate arrays [f,v]. It is not critical, just a choice to ease the next call to the patch object.

I need to retrieve the vertices coordinates:

%// Retrieve the vertices coordinates
X = fv.vertices(:,1) ;
Y = fv.vertices(:,2) ;
Z = fv.vertices(:,3) ;

You can then apply your transformation. I choose a simple one in this example, but any transformation function is valid.

%// Transform
X = -X.*Y.^2 ;
Y = Y.*X ;
Z = Z*2 ;

Then I rebuild a new structure for the patch which will display the transformed object.
This is the important bit:

%// create new patch structure
fvt.vertices = [X Y Z] ;   %// with the new transformed 'vertices'
fvt.faces = fv.faces ;     %// but we keep the same 'faces'

Then I display it the same way (well with a slightly different angle for a better view):

%// Plot the transformed isosurface
figure(2);cla
pt = patch( fvt ,'FaceColor','red','EdgeColor','none');
%// refine the view
grid off ; set(gca,'Color','none') ; daspect([1,1,1]) ; view(-3,4) ; axis tight ; camlight ; lighting gouraud

Which produces the figure:
Transformed isosurface

(If you paste all the code snippet in one file it should run and give you the same output.)

Hoki
  • 11,637
  • 1
  • 24
  • 43