Return to Answer

Rather than fit bilinear patches to your sample points, which are only C0 continuous$^\spadesuit$ may I instead suggest fitting the surface using Catmull-Rom splines? These will give you C2 continuity. Catmull-Rom is an interpolating spline which means it will pass through your sample points, as did the bilinear surface, but will do so smoothly from patch to patch.

$^\spadesuit$ because the first derivatives won't generally match at the boundaries of your patches

To form each bicubic Catmull-Rom patch, $C_{i,j}$, which has a given 2x2 set of samples $\{S_{i,j}, S_{i+1,j}, S_{i,j+1}, S_{i+1,j+1}\}$ at the corners, you need a set of 4x4 samples formed from the 2x2 at the centre and the surrounding 12 samples, i.e.$\{S_{i-1,j-1}$ through $S_{i+2,j+2}\}$, as shown below. (My apologies that this is 'hand placed' so is far from ideal) [![enter image description here][1]][1]

If we just consider the behaviour of the surface at sample, $S_{i,j}$, then the partial derivative in the $i$ direction is given by $$\frac{\partial C_{ij}}{\partial i}=\frac{S_{i+1,j}-S_{i-1,j}}{2}$$ and similarly, in the $j$ direction by $$\frac{\partial C_{ij}}{\partial j}=\frac{S_{i,j+1}-S_{i,j-1}}{2}$$

Remembering these are vectors, you can use the cross product to generate a surface normal at any of the $S_{i,j}$ points$^\clubsuit$. (Though you may need an extra set of boundary samples)

Ali: You mentioned that both partial derivatives are vectors(in i,j directions). could you write each vector components plz

Certainly: If we say $S_{i,j}=[x_{i,j},y_{i,j},z_{i,j}]^T$ then $$\frac{\partial C_{ij}}{\partial i}=\frac{1}{2}.[(x_{i+1,j}-x_{i-1,j}),(y_{i+1,j}-y_{i-1,j}),(z_{i+1,j}-z_{i-1,j})]^T$$ $$\frac{\partial C_{ij}}{\partial j}=\frac{1}{2}.[(x_{i,j+1}-x_{i,j-1}),(y_{i,j+1}-y_{i,j-1}),(z_{i,j+1}-z_{i,j-1})]^T$$

<br>If you've sampled the surface 'finely' enough you could probably choose to render the surface by just using 2 or 4 triangles per original patch.

Alternatively, you can recursively subdivide your quad mesh to get a finer set of triangles. To double the number of points in the 'i' direction, simply insert a new point between each $S_{i,j}$ and $S_{i+1,j}$ that equals ... $$\frac{9(S_{i,j}+S_{i+1,j}) - (S_{i-1,j}+S_{i+2,j})}{16}$$ Once you have the finer grid in 'i' simply repeat with 'j' with the appropriate equivalent set of points. You can repeat this ad nauseam to get whatever triangle resolution you desire ;-)

$^\clubsuit$ (Copying some text from the comments below): I've made a couple of valid assumptions here. Noting that if a given first (partial) derivative is non-zero (which is guaranteed in your case because of the way you construct the surface), then it also is locally tangential to that surface. Further, we know your two partial derivates aren't parallel, so the cross product will construct a valid normal. FWIW For some arrangements of control points, i.e. in more general use cases, you can get so-called 'degenerate' surfaces where you can't naively apply the above rules, but you're safe here. [1]: https://i.sstatic.net/M48A3.png

Rather than fit bilinear patches to your sample points, which are only C0 continuous$^\spadesuit$ may I instead suggest fitting the surface using Catmull-Rom splines? These will give you C2 continuity. Catmull-Rom is an interpolating spline which means it will pass through your sample points, as did the bilinear surface, but will do so smoothly from patch to patch.

$^\spadesuit$ because the first derivatives won't generally match at the boundaries of your patches

Background/Assumption

You have a function $$F:\Re^2\to\Re^3$$ of the form $$F(s,t)=[s,t,f(s,t)]$$ where $f:\Re^2\to\Re$. The function, $f$, can be considered a 'black box' in the sense that we might not be given it explicitly.

$F$ has been regularly sampled as $$S_{i,j}=F(x_{scale}.i, y_{scale}.j, f(x_{scale}.i,y_{scale}.j))$$ where $i,j\in\mathbb{Z}$. You then want interpolate between these sample points using piecewise bicubic Catmull-Rom patches.

Constructing patches

To form a particular bicubic Catmull-Rom patch, $C_{i,j}$, which has a given 2x2 set of samples $\{S_{i,j}, S_{i+1,j}, S_{i,j+1}, S_{i+1,j+1}\}$ at the corners, you need the set of 4x4 samples formed from the 2x2 at the centre and the surrounding 12 samples, i.e.$\{S_{i-1,j-1}$ through $S_{i+2,j+2}\}$, as shown below. (My apologies that this is 'hand placed' so is far from ideal) [![enter image description here][1]][1]

If we just consider the behaviour of the surface at sample, $S_{i,j}$, then the partial derivative in the $i$ direction is given by $$\frac{\partial C_{ij}}{\partial i}=\frac{S_{i+1,j}-S_{i-1,j}}{2}$$ and similarly, in the $j$ direction by $$\frac{\partial C_{ij}}{\partial j}=\frac{S_{i,j+1}-S_{i,j-1}}{2}$$

Remembering these are vectors, you can use the cross product to generate a surface normal at any of the $S_{i,j}$ points$^\clubsuit$. (Though you may need an extra set of boundary samples)

Ali: You mentioned that both partial derivatives are vectors(in i,j directions). could you write each vector components plz

Certainly: If we say $S_{i,j}=[x_{i,j},y_{i,j},z_{i,j}]^T$ then $$\frac{\partial C_{ij}}{\partial i}=\frac{1}{2}.[(x_{i+1,j}-x_{i-1,j}),(y_{i+1,j}-y_{i-1,j}),(z_{i+1,j}-z_{i-1,j})]^T$$ $$\frac{\partial C_{ij}}{\partial j}=\frac{1}{2}.[(x_{i,j+1}-x_{i,j-1}),(y_{i,j+1}-y_{i,j-1}),(z_{i,j+1}-z_{i,j-1})]^T$$

<br>If you've sampled the surface 'finely' enough you could probably choose to render the surface by just using 2 or 4 triangles per original patch.

Alternatively, you can recursively subdivide your quad mesh to get a finer set of triangles. To double the number of points in the 'i' direction, simply insert a new point between each $S_{i,j}$ and $S_{i+1,j}$ that equals ... $$\frac{9(S_{i,j}+S_{i+1,j}) - (S_{i-1,j}+S_{i+2,j})}{16}$$ Once you have the finer grid in 'i' simply repeat with 'j' with the appropriate equivalent set of points. You can repeat this ad nauseam to get whatever triangle resolution you desire ;-)

$^\clubsuit$ (Copying some text from the comments below): I've made a couple of valid assumptions here. Noting that if a given first (partial) derivative is non-zero (which is guaranteed in your case because of the way you construct the surface), then it also is locally tangential to that surface. Further, we know your two partial derivates aren't parallel, so the cross product will construct a valid normal. FWIW For some arrangements of control points, i.e. in more general use cases, you can get so-called 'degenerate' surfaces where you can't naively apply the above rules, but you're safe here. [1]: https://i.sstatic.net/M48A3.png

Rather than fit bilinear patches to your sample points, which are only C0 continuous$^\spadesuit$ may I instead suggest fitting the surface using Catmull-Rom splines? These will give you C2 continuity. Catmull-Rom is an interpolating spline which means it will pass through your sample points, as did the bilinear surface, but will do so smoothly from patch to patch.

$^\spadesuit$ because the first derivatives won't generally match at the boundaries of your patches

To form eachbicubic Catmull-Rom patch, $C$, which has a given 2x2 set of samples $\{S_{i,j}, S_{i+1,j}, S_{i,j+1}, S_{i+1,j+1}\}$ at the corners, you need a set of 4x4 samples formed from the 2x2 at the centre and the surrounding 12 samples, i.e.$\{S_{i-1,j-1}$ through $S_{i+2,j+2}\}$, as shown below. (My apologies that this is 'hand placed' so is far from ideal)

If we just consider the sample, $S_{i,j}$, then the partial derivative in the 'i' direction is given by $$\frac{\partial C}{\partial i}=\frac{S_{i+1,j}-S_{i-1,j}}{2}$$ and similarly, in the 'j' direction by $$\frac{\partial C}{\partial j}=\frac{S_{i,j+1}-S_{i,j-1}}{2}$$

Remembering these are vectors, you can use the cross product to generate a surface normal at any of the $S_{i,j}$ points. (Though you may need an extra set of boundary samples)

Ali: You mentioned that both partial derivatives are vectors(in i,j directions). could you write each vector components plz

Certainly: If we say $S_{i,j}=[x_{i,j},y_{i,j},z_{i,j}]^T$ then $$\frac{\partial C}{\partial i}=\frac{1}{2}.[(x_{i+1,j}-x_{i-1,j}),(y_{i+1,j}-y_{i-1,j}),(z_{i+1,j}-z_{i-1,j})]^T$$ $$\frac{\partial C}{\partial j}=\frac{1}{2}.[(x_{i,j+1}-x_{i,j-1}),(y_{i,j+1}-y_{i,j-1}),(z_{i,j+1}-z_{i,j-1})]^T$$

<br>If you've sampled the surface 'finely' enough you could probably choose to render the surface by just using 2 or 4 triangles per original patch.

Alternatively, you can recursively subdivide your quad mesh to get a finer set of triangles. To double the number of points in the 'i' direction, simply insert a new point between each $S_{i,j}$ and $S_{i+1,j}$ that equals ... $$\frac{9(S_{i,j}+S_{i+1,j}) - (S_{i-1,j}+S_{i+2,j})}{16}$$ Once you have the finer grid in 'i' simply repeat with 'j' with the appropriate equivalent set of points. You can repeat this ad nauseam to get whatever triangle resolution you desire ;-)

Rather than fit bilinear patches to your sample points, which are only C0 continuous$^\spadesuit$ may I instead suggest fitting the surface using Catmull-Rom splines? These will give you C2 continuity. Catmull-Rom is an interpolating spline which means it will pass through your sample points, as did the bilinear surface, but will do so smoothly from patch to patch.

$^\spadesuit$ because the first derivatives won't generally match at the boundaries of your patches

To form each bicubic Catmull-Rom patch, $C_{i,j}$, which has a given 2x2 set of samples $\{S_{i,j}, S_{i+1,j}, S_{i,j+1}, S_{i+1,j+1}\}$ at the corners, you need a set of 4x4 samples formed from the 2x2 at the centre and the surrounding 12 samples, i.e.$\{S_{i-1,j-1}$ through $S_{i+2,j+2}\}$, as shown below. (My apologies that this is 'hand placed' so is far from ideal) [![enter image description here][1]][1]

If we just consider the behaviour of the surface at sample, $S_{i,j}$, then the partial derivative in the $i$ direction is given by $$\frac{\partial C_{ij}}{\partial i}=\frac{S_{i+1,j}-S_{i-1,j}}{2}$$ and similarly, in the $j$ direction by $$\frac{\partial C_{ij}}{\partial j}=\frac{S_{i,j+1}-S_{i,j-1}}{2}$$

Remembering these are vectors, you can use the cross product to generate a surface normal at any of the $S_{i,j}$ points$^\clubsuit$. (Though you may need an extra set of boundary samples)

Ali: You mentioned that both partial derivatives are vectors(in i,j directions). could you write each vector components plz

Certainly: If we say $S_{i,j}=[x_{i,j},y_{i,j},z_{i,j}]^T$ then $$\frac{\partial C_{ij}}{\partial i}=\frac{1}{2}.[(x_{i+1,j}-x_{i-1,j}),(y_{i+1,j}-y_{i-1,j}),(z_{i+1,j}-z_{i-1,j})]^T$$ $$\frac{\partial C_{ij}}{\partial j}=\frac{1}{2}.[(x_{i,j+1}-x_{i,j-1}),(y_{i,j+1}-y_{i,j-1}),(z_{i,j+1}-z_{i,j-1})]^T$$

<br>If you've sampled the surface 'finely' enough you could probably choose to render the surface by just using 2 or 4 triangles per original patch.

Alternatively, you can recursively subdivide your quad mesh to get a finer set of triangles. To double the number of points in the 'i' direction, simply insert a new point between each $S_{i,j}$ and $S_{i+1,j}$ that equals ... $$\frac{9(S_{i,j}+S_{i+1,j}) - (S_{i-1,j}+S_{i+2,j})}{16}$$ Once you have the finer grid in 'i' simply repeat with 'j' with the appropriate equivalent set of points. You can repeat this ad nauseam to get whatever triangle resolution you desire ;-)

$^\clubsuit$ (Copying some text from the comments below): I've made a couple of valid assumptions here. Noting that if a given first (partial) derivative is non-zero (which is guaranteed in your case because of the way you construct the surface), then it also is locally tangential to that surface. Further, we know your two partial derivates aren't parallel, so the cross product will construct a valid normal. FWIW For some arrangements of control points, i.e. in more general use cases, you can get so-called 'degenerate' surfaces where you can't naively apply the above rules, but you're safe here. [1]: https://i.sstatic.net/M48A3.png

Rather than fit bilinear patches to your sample points, which are only C0 continuous$^\spadesuit$ may I instead suggest fitting the surface using Catmull-Rom splines, which will give you C2 continuity? Catmull-Rom is an interpolating spline which means it will pass through your sample points as did the bilinear surface, but will do so smoothly.

$^\spadesuit$ because the first derivatives won't generally match at the boundaries of your patches

To form a bicubic Catmull-Rom patch which has a given 2x2 set of samples $\{S_{i,j}, S_{i+1,j}, S_{i,j+1}, S_{i+1,j+1}\}$ at the corners, you need a set of 4x4 samples formed from the 2x2 at the centre and the surrounding 12 samples, i.e.$\{S_{i-1,j-1}$ through $S_{i+2,j+2}\}$, as shown below. (My apologies that this is 'hand placed' so is far from ideal)

If we just consider the sample, $S_{i,j}$, then the partial derivative in the 'i' direction is given by $$\frac{S_{i+1,j}-S_{i-1,j}}{2}$$ and similarly, in the 'j' direction by $$\frac{S_{i,j+1}-S_{i,j-1}}{2}$$

Remembering these are vectors, you can use the cross product to generate a surface normal at any of the $S_{i,j}$ points. (Though you may need an extra set of boundary samples)

If you've sampled the surface 'finely' enough you could probably choose to render the surface by just using 2 or 4 triangles per original patch.

Alternatively, you can recursively subdivide your quad mesh to get a finer set of triangles. To double the number of points in the 'i' direction, simply insert a new point between each $S_{i,j}$ and $S_{i+1,j}$ that equals ... $$\frac{9(S_{i,j}+S_{i+1,j}) - (S_{i-1,j}+S_{i+2,j})}{16}$$ Once you have the finer grid in 'i' simply repeat with 'j' with the appropriate equivalent set of points. You can repeat this ad nauseam to get whatever triangle resolution you desire ;-)

Rather than fit bilinear patches to your sample points, which are only C0 continuous$^\spadesuit$ may I instead suggest fitting the surface using Catmull-Rom splines? These will give you C2 continuity. Catmull-Rom is an interpolating spline which means it will pass through your sample points, as did the bilinear surface, but will do so smoothly from patch to patch.

$^\spadesuit$ because the first derivatives won't generally match at the boundaries of your patches

To form eachbicubic Catmull-Rom patch, $C$, which has a given 2x2 set of samples $\{S_{i,j}, S_{i+1,j}, S_{i,j+1}, S_{i+1,j+1}\}$ at the corners, you need a set of 4x4 samples formed from the 2x2 at the centre and the surrounding 12 samples, i.e.$\{S_{i-1,j-1}$ through $S_{i+2,j+2}\}$, as shown below. (My apologies that this is 'hand placed' so is far from ideal)

If we just consider the sample, $S_{i,j}$, then the partial derivative in the 'i' direction is given by $$\frac{\partial C}{\partial i}=\frac{S_{i+1,j}-S_{i-1,j}}{2}$$ and similarly, in the 'j' direction by $$\frac{\partial C}{\partial j}=\frac{S_{i,j+1}-S_{i,j-1}}{2}$$

Remembering these are vectors, you can use the cross product to generate a surface normal at any of the $S_{i,j}$ points. (Though you may need an extra set of boundary samples)

Ali: You mentioned that both partial derivatives are vectors(in i,j directions). could you write each vector components plz

Certainly: If we say $S_{i,j}=[x_{i,j},y_{i,j},z_{i,j}]^T$ then $$\frac{\partial C}{\partial i}=\frac{1}{2}.[(x_{i+1,j}-x_{i-1,j}),(y_{i+1,j}-y_{i-1,j}),(z_{i+1,j}-z_{i-1,j})]^T$$ $$\frac{\partial C}{\partial j}=\frac{1}{2}.[(x_{i,j+1}-x_{i,j-1}),(y_{i,j+1}-y_{i,j-1}),(z_{i,j+1}-z_{i,j-1})]^T$$

<br>If you've sampled the surface 'finely' enough you could probably choose to render the surface by just using 2 or 4 triangles per original patch.

Alternatively, you can recursively subdivide your quad mesh to get a finer set of triangles. To double the number of points in the 'i' direction, simply insert a new point between each $S_{i,j}$ and $S_{i+1,j}$ that equals ... $$\frac{9(S_{i,j}+S_{i+1,j}) - (S_{i-1,j}+S_{i+2,j})}{16}$$ Once you have the finer grid in 'i' simply repeat with 'j' with the appropriate equivalent set of points. You can repeat this ad nauseam to get whatever triangle resolution you desire ;-)