Saving limit surfaces

[Gregory Tada]

Hey guys,

I haven't programmed anything for 20 years, but I want the ability in Cinema 4D to save out limit surfaces as an IGES file, so I started creating a plug-in. I'm using the libIGES library to handle the file output. 

I'm looping through the patchtable and saving out each 4x4 surface (I think). The overall shape of the result looks like I'm going in the right direction, but the patches don't seem to have the correct control points. 

I have attached the IGES file and a screenshot of it if anybody would care to look at it. 

If anybody can point me in the right direction I would deeply appreciate it! Thank you in advance! 

Cheers!

Greg

Here's a screenshot of what my code produces with the example pyramid mesh: 

Here is my code: 

#define _USE_MATH_DEFINES // for C++

#include <cmath>

#include <cassert>

#include <cstring>

#include <cfloat>

#include <api/dll_iges.h>

#include <api/all_api_entities.h>

#include <opensubdiv/far/topologyDescriptor.h>

#include <opensubdiv/far/primvarRefiner.h>

#include <opensubdiv/far/patchTableFactory.h>

#include <opensubdiv/far/patchMap.h>

#include <opensubdiv/far/ptexIndices.h>

#include "helpers.h"

using namespace OpenSubdiv;

static int const g_nverts = 5;

static double const g_verts[24] = { 0.0f,   0.0f, 20.0f,

                                    0.0f, -20.0f,  0.0f,

                                   20.0f,   0.0f,  0.0f,

                                    0.0f,  20.0f,  0.0f,

                                  -20.0f,   0.0f,  0.0f, };

static int const g_vertsperface[5] = { 3, 3, 3, 3, 4 };

static int const g_nfaces = 5;

static int const g_faceverts[16] = { 0, 1, 2,

                                     0, 2, 3,

                                     0, 3, 4,

                                     0, 4, 1,

                                     4, 3, 2, 1 };

static int const g_ncreases = 4;

static int const g_creaseverts[8] = { 4, 3, 3, 2, 2, 1, 1, 4 };

static float const g_creaseweights[4] = { 3.0f, 3.0f, 3.0f, 3.0f };

// Creates a Far::TopologyRefiner from the pyramid shape above

static Far::TopologyRefiner* createTopologyRefiner();

//------------------------------------------------------------------------------

int main(int, char**) {

    // Generate a FarTopologyRefiner (see far_tutorial_0 for details).

    Far::TopologyRefiner* refiner = createTopologyRefiner();

    // Adaptively refine the topology with an isolation level capped at 3

    // because the sharpest crease in the shape is 3.0f (in g_creaseweights[])

    int maxIsolation = 3;

    refiner->RefineAdaptive(

        Far::TopologyRefiner::AdaptiveOptions(maxIsolation));

    // Generate a set of Far::PatchTable that we will use to evaluate the

    // surface limit

    Far::PatchTableFactory::Options patchOptions;

    patchOptions.endCapType =

        Far::PatchTableFactory::Options::ENDCAP_GREGORY_BASIS;

    Far::PatchTable const* patchTable =

        Far::PatchTableFactory::Create(*refiner, patchOptions);

    // Compute the total number of points we need to evaluate patchtable.

    // we use local points around extraordinary features.

    int nRefinerVertices = refiner->GetNumVerticesTotal();

    int nLocalPoints = patchTable->GetNumLocalPoints();

    // Create a buffer to hold the position of the refined verts and

    // local points, then copy the coarse positions at the beginning.

    std::vector<Vertex> verts(nRefinerVertices + nLocalPoints);

    memcpy(&verts[0], g_verts, g_nverts * 3 * sizeof(double));

    // Adaptive refinement may result in fewer levels than maxIsolation.

    int nRefinedLevels = refiner->GetNumLevels();

    // Interpolate vertex primvar data : they are the control vertices

    // of the limit patches (see far_tutorial_0 for details)

    Vertex* src = &verts[0];

    for (int level = 1; level < nRefinedLevels; ++level) 

    {

        Vertex* dst = src + refiner->GetLevel(level - 1).GetNumVertices();

        Far::PrimvarRefiner(*refiner).Interpolate(level, src, dst);

        src = dst;

    }

    // Evaluate local points from interpolated vertex primvars.

    patchTable->ComputeLocalPointValues(&verts[0], &verts[nRefinerVertices]);

    // instantiate the IGES data object ///////////////////////////////////////

    DLL_IGES model;

    // Create knot vector

    double knots[] = { 0.0, 0.0, 0.0, 0.0, 1.0, 1.0, 1.0, 1.0 };

    // Create a surface object

    DLL_IGES_ENTITY_128 surf(model, true);

    // IGES object created ////////////////////////////////////////////////////

    // Loop through each patch and save out 4x4 vertices each

    int na = patchTable->GetNumPatchArrays();

    bool error;

    for (int i = 0; i < na; i++)

    {

        Far::PatchDescriptor pd = patchTable->GetPatchArrayDescriptor(i);

        if (pd == 6) // Type::REGULAR

        {

            Far::ConstIndexArray arraycvs = patchTable->GetPatchArrayVertices(i);

            int np = patchTable->GetNumPatches(i);

            for (int patch = 0; patch < np; patch++)

            {

                Far::ConstIndexArray cvs = patchTable->GetPatchVertices(i, patch);

                int cvCount = cvs.size();

                std::vector<double> surfVerts;

                for (int cv = 0; cv < cvCount; cv++)

                {

                    surfVerts.push_back(verts[cvs[cv]].point[0]);

                    surfVerts.push_back(verts[cvs[cv]].point[1]);

                    surfVerts.push_back(verts[cvs[cv]].point[2]);

                }

                if (patch != 0)

                    error = surf.NewEntity();

                error = surf.SetNURBSData(4, 4, 4, 4, knots, knots, &surfVerts[0], false, false, false, 0.0, 1.0, 0.0, 1.0);

                surfVerts.clear();

            }

        }

    }

    model.Write("test.igs", true);

    return 0;

}

static Far::TopologyRefiner* createTopologyRefiner() 

{

    typedef Far::TopologyDescriptor Descriptor;

    Sdc::SchemeType type = OpenSubdiv::Sdc::SCHEME_CATMARK;

    Sdc::Options options;

    options.SetVtxBoundaryInterpolation(Sdc::Options::VTX_BOUNDARY_EDGE_ONLY);

    Descriptor desc;

    desc.numVertices = g_nverts;

    desc.numFaces = g_nfaces;

    desc.numVertsPerFace = g_vertsperface;

    desc.vertIndicesPerFace = g_faceverts;

    desc.numCreases = g_ncreases;

    desc.creaseVertexIndexPairs = g_creaseverts;

    desc.creaseWeights = g_creaseweights;

    // Instantiate a FarTopologyRefiner from the descriptor.

    Far::TopologyRefiner* refiner =

        Far::TopologyRefinerFactory<Descriptor>::Create(desc,

            Far::TopologyRefinerFactory<Descriptor>::Options(type, options));

    return refiner;

}

using namespace std;

[Gregory Tada]

Ah, figured out that my knot vector was the issue: 

{ 0.0, 0.0, 0.0, 0.0, 1.0, 1.0, 1.0, 1.0 };

What I needed was: 

{ 0.0, 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0 };

Cheers!

Greg

[David G Yu]

 Nice! Good to hear that you got it working!

-David

On Sep 27, 2021, at 7:48 PM, Gregory Tada <greg...@gmail.com> wrote:

Ah, figured out that my knot vector was the issue: 

--
You received this message because you are subscribed to the Google Groups "OpenSubdiv Forum" group.
To unsubscribe from this group and stop receiving emails from it, send an email to opensubdiv+...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/opensubdiv/32907ebf-4560-4045-b100-cc3a60ef4063n%40googlegroups.com.

[Gregory Tada]

Thanks David! However, it's not quite working 100%, and I'm wondering if anybody can point me in the right direction. The Gregory Basis patches... I don't understand them. They're order 4 in one direction and order 5 in the other? 20 control points? How do I calculate these patches? I've tried a few experiments and gotten some strange results (see attached). Any guidance is welcome and appreciated! 

Greg

[David G Yu]

You can think of a Gregory Basis patch as an extension of a bicubic Bezier patch.

The additional interior control points allow the patch to satisfy tangent plane continuity along the boundary edges of the patch. The organization of the control points may be surprising initially, but follows from the sequence of points at the 4 corners of the patch. See the diagram below.

 

Within the OpenSubdiv code base, evaluation of Gregory Basis patches is implemented in patchBasis.cpp#L348 as described in Loop, Schaefer, Ni, Castano SIGGRAPH Asia 2009

-David

To view this discussion on the web visit https://groups.google.com/d/msgid/opensubdiv/9d54bd9e-d88c-40d6-a3fb-5522c00d3094n%40googlegroups.com.

[Gregory Tada]

David,

That was extremely helpful, thank you! I cheated by picking just one of the paired control points (4, 9, 14, 19). For my purposes, I wonder if averaging the pairs would be good enough (i.e., (f0+ + f0-) / 2). It's a hack, but it's enough to fill that surface and move on to creating the plugin code for Cinema 4D. 

Cheers!

Greg

[Pieter Barendrecht]

Hi Greg,

The rationale of using Gregory patches is to fill the regions in the (limit) surface around the so-called extraordinary vertices (i.e. vertices with a valency — that is, the number of edges connecting to that vertex — different from 4) in a smooth way. Without getting too technical (hopefully), this means that Gregory patches connect to their neighbouring (Bézier/B-spline) patches not only continuously (i.e. without gaps), but actually in a way such that the tangent plane varies continuously. This latter part means that you won't get 'kinks' in your surface (in a curve setting, an example of a kink would be the behaviour of the absolute function at the origin — it is continuous, but the line indicating the slope of the function (i.e. the tangent line) changes rather suddenly at the origin, i.e. it does not change continuously). Visually, these kinks are typically something you'd like to avoid (as they tend to be rather noticeable, especially if the surface reflects its environment).

Now, if you pick either just one vertex for each pair of inner control points (i.e. your P4, P9, P14 and P19) or use some pairwise average, the resulting patches will merely connect continuously to their neighbours. That is, they will generally not connect smoothly, which means that you might get these kinks in your surface (though in some cases, e.g. for meshes with local symmetries, they might not be noticeable).

Looking at the Gregory patch in more detail, the pairs of inner control points can be explained by looking at (bicubic) Bézier patches. For two such patches — let's call them L and R — to connect smoothly, the pairs of interior points to either side of the boundary curve shared by L and R (think : | : with | the shared boundary and : a pair of interior points) have to satisfy certain conditions. No problem thus far. However, if you now want to smoothly attach a third Bézier patch B, say below L, there is another set of conditions that should be satisfied, which now involves the pairs of interior points to either side of the boundary curve shared by L and B:

L  R

B

But note that one of the interior points of patch L involved here — namely the bottom-right one — is also involved in the smoothness conditions between L and R... Ultimately, it turns out that in general, the four interior points of the (bicubic) Bézier patch are not sufficient to guarantee smoothly connecting neighbours. The Gregory patch solves this by using two interior points per (shared) boundary edge, which leads to those eight points from David's figure.

Finally, the boundary of a Gregory patch is identical to that of a (bicubic) Bézier patch (that is, the inner points are not involved). However, in the interior, there is a clear difference — the Gregory patch uses linear (or well, affine) combinations of the pairs of inner control points that depend on the (u,v) values. For example, for the pair (P3, P4), we get (uP3 + vP4)/(u + v). In other words, for each (u,v) value, you will obtain four interior points, each a linear (affine) combination of one of those pairs. From that point of view, the patch 'locally' acts as a (bicubic) Bézier patch.

Cheers,

Pieter

To view this discussion on the web visit https://groups.google.com/d/msgid/opensubdiv/696bd665-2eec-4f1a-8d25-93ae93bac52cn%40googlegroups.com.

[Gregory Tada]

Pieter, 

You kinda broke my brain there. ;) I really appreciate the explanation; I may have to digest it for a while to fully understand. So those inner points act almost like levers on the associated edge and only that edge? 

The challenge for me is that no NURBS modeler that I use supports Gregory patches (Autodesk Alias and McNeel Rhino), so if I want to export a subdivision surface model from Cinema 4D into Alias for example, I have to find a reasonably accurate and quick way to convert those Gregory patches into B-spline patches. I've heard that one of Autodesk's apps approximates a Gregory patch with a 7th degree/order 8 surface, but I'd like to keep everything 5th degree/order 6 or lower. 

I did open up the resulting surfaces in Alias, and they definitely break continuity: 

In the end I'll probably have to resort to approximating the Gregory patches. For my intended use, that might be good enough. Thank you both! 

Cheers!

Greg

[David G Yu]

If approximating the surface and tolerating tangent discontinuities is acceptable for your use case then you can have OpenSubdiv generate an approximation for you via the PatchTableFactory::Options.

i.e. if you specify:

    patchOptions.endCapType = Far::PatchTableFactory::Options::ENDCAP_BSPLINE_BASIS;

instead of:

    patchOptions.endCapType = Far::PatchTableFactory::Options::ENDCAP_GREGORY_BASIS;

in your example code, then irregular regions of the surface will be approximated using BSpline patches instead of Gregory Basis patches.

-David

To view this discussion on the web visit https://groups.google.com/d/msgid/opensubdiv/21f07ee2-ba79-4452-a059-ee1b7ad6341cn%40googlegroups.com.

[Pieter Barendrecht]

Hi Greg,

This goes a bit outside the scope of OpenSubdiv, but I wanted to mention it nevertheless — Gregory patches can only be converted to rational patches, not to polynomial ones. More specifically, the type of Gregory patch mentioned above (there are also other types of Gregory patches) can be converted to a degree 7x7 rational Bézier patch, which can then also be represented as a degree 7x7 NURBS patch (the R in NURBS stands for rational). Exact conversion of Gregory patches to polynomial Bézier- or B-spline patches is not possible.

However, seeing as you would like to work with patches that are at most bi-degree 5 (= bi-quintic), you will have to look for an alternative construction (not based on Gregory patches) that results in patches connecting smoothly to each other (keywords here are 'tangent continuous', also known as 'G1 continuity'). If I remember correctly, the modern T-spline surfaces use a polynomial bi-degree 4 or 5 construction to smoothly fill the regions around extraordinary vertices (referred to as 'star points' in T-spline terminology). I'll be working on a related problem for an upcoming research project, so I'd be happy to help out with some details if time allows it. Just drop me a line :).

Cheers,

Pieter

To view this discussion on the web visit https://groups.google.com/d/msgid/opensubdiv/CACszWpZHJ25L06w3ehg-%3DNr2v2PknRhjonfQ%2B83su3epu%2Bq1UQ%40mail.gmail.com.

[Gregory Tada]

David,

Thank you for the suggestion; I may do that. For now I'll backburner the Gregory patches and get the rest of the code going until I figure out a better solution for my use case. 

Cheers!

Greg

[Gregory Tada]

Pieter, 

I might go for the 7x7 Bézier conversion in the future if my current solution is deemed unsuitable. I've found a few online articles on the topic, but if you know of an easier-to-digest version, please let me know. :D 

I would definitely like to know more about your upcoming research project. I'll contact you about it. 

Thank you guys for the help! I'm a rank novice and really appreciate the guidance. 

Cheers,

Greg

[Gregory Tada]

I was able to get the Cinema 4D plug-in portion to work, and it now converts a base mesh and saves that to IGES. Thank you guys for the help! 

One issue I'm running into: when I give it an open mesh, the resulting surfaces are odd. I took a cube and deleted one face: 

I'm going to dig in and double-check my code, but anybody have an idea where to look? 

Thanks!

Greg

[David G Yu]

The patches along sharp boundaries (and along inf sharp creases, which are a kind of interior boundary edge) require special treatment.

The OpenSubdiv PatchTable defines bits in PatchParam to allow sharp boundaries of parametric patches to be evaluated correctly. See far/patchParam.h and far/patchBasis.cpp along with diagrams at patch-parameterization

The calculation of the necessary "phantom" control points is pretty straightforward (see the links above), and I think you should be able to account for them in your export plugin.

Thanks!

-David

To view this discussion on the web visit https://groups.google.com/d/msgid/opensubdiv/bbf2a368-fa72-4cdd-a63f-c9b2a442b021n%40googlegroups.com.

[Gregory Tada]

Got it working. :D Thanks go out to David and Pieter for the excellent guidance! 

[David G Yu]

Looking good!

-David

To view this discussion on the web visit https://groups.google.com/d/msgid/opensubdiv/ee43a347-420f-4404-a5ab-2d4c4df1679an%40googlegroups.com.

[Nick Kallen]

Hi - I know this is an old thread but I was wondering if anyone could link to resources on converting gregory patches to 7x7 rational b-spline patches?

Thanks!

[Pieter Barendrecht]

Hi Nick,

Here's a link to a paper discussing the conversion to rational Bézier patches (which you could then convert to B-spline patches) — https://doi.org/10.1007/978-4-431-68123-6_32. Farin also briefly describes it in his book 'NURBS for Curve and Surface Design'. Hope that helps!

Best,

Pieter

To view this discussion on the web visit https://groups.google.com/d/msgid/opensubdiv/1688eae9-6fc1-472b-99ba-8af44e760fden%40googlegroups.com.