Hello,

I'm a trainee for software development and in my free time I try to do various things. I decided I wanted to figure out how "Dual Contouring" works, but I haven't been successful yet. I could find some implementations in code which are very hard to follow along and same sheets of paper "explaining" it written in a way that I have a hard time even understanding what the input and output is. All I know that it is used to make  a surface/mesh out of voxels. Is there any explanation someone can understand without having studied math/computer science? The problem is also that most of the words I can't even translate to German(like Hermite Data) nor I can find a explanation which I can understand. I know how Marching Cubes work but I just can't find anything to Dual Contouring also I don't quite get the sense of octrees. As far I'm aware of this is a set of data organized like a tree where each element points to 8 more elements. I don't get how that could be helpful in a infinite 3D world and why I shouldn't just use a List of coordinates with zeros and ones for dirt/air or something like that.

Thanks

A clueless trainee ^^

Have you seen: https://github.com/emilk/Dual-Contouring ?

It compiles and runs just fine for me on MacOS X. No such luck compiling when using Visual Studio (cl).  

I too have been looking forward to familiarizing myself with Dual Contouring.

I have looked through the code, but I have some issues with the variable names and functions used. As far what I have seen the Field consists of a array of planes, what is kinda confusing as I thought Dual Contouring is only using like a 3D array of zeros and ones. Also I'm not quite sure if a voxel is in the middle of a cube or on the corner because it seems like it is iterating with the foreach3D(function?) through each voxel and determines if the surrounding voxels are all outside or inside(below zero and greater zero as it seems). If thats not the case then it makes vertexes to CrossingCorners as it seems. Is that correct for the first part? I find this confusing because voxels in Minecraft terms are simply displaying 1 Cube and can either be 0 or a ID of a block, but here it seems  to be that 8 voxels/coordinates represent 1 Cube or whatever might come out of the algorithm.

Have you tried emailing the author of that code? emil.ernerfeldt@gmail.com

39 minutes ago, sjhalayka said:

Have you tried emailing the author of that code? emil.ernerfeldt@gmail.com

Just did so, I hope he is willing to help me ^^

Cool. Let me know how it goes.

25 minutes ago, swiftcoder said:

Conceptually, Dual Contouring is not that different from Marching Cubes. In both cases, you divide your volume data into a grid, and at each cube in the grid you evaluate the surface. Any cube containing edge-crossing - or in other words, any cube that intersects the surface (i.e. isn't entirely empty or entirely full) will generate output geometry.

The key difference between the two is that in Marching Cubes you will just directly generate geometry from the edge crossings. In Dual Contouring you will instead use the normals at each edge-crossing (this is why you need Hermite Data) to locate a single exact point on the surface, and you'll generate geometry from these points (spanning multiple cubes) in a subsequent pass.

The complicated part of Dual Contouring is mostly in how you find the optimal point based on the normals of neighbouring edge crossings. The original paper uses QEF minimisation, which is... a little opaque to those of us without a deep math background.

If you are looking for a more readable implementation, Lin has a pretty simple one on github.

Wow, thank you so much you just cleared up so much stuff I found to be confusing. I don't quite understand how those normals should be working. If for example in my game a Cube gets placed by the player then I will create volume data containing the coordinates of the 8 edges, right? But what normal direction does I give the 8 coordinates and how exactly do I use it? Also how should I locate a point on the surface? I thought Dual Contouring is creating this surface out of a set of data, so how can I determine a location on a shape then?

I made a image on how I thought it is working:

1. User wants to create a cube

2. create 8 voxels

3. Add/Calculate Normals(free hand ^^)

4.Dual Contouring 

Is this how it works or did I miss something? And how does Dual Contouring actually use those normals?

I took a shot at deciphering the emilk code, and I'm sort of half-triumphant.

I put up a copy of the code in my own GitHub repository:

https://github.com/sjhalayka/dc

I added a function in main.cpp called addFractal (also see below for source code), which adds a quaternion Julia set. I set the quaternion parameters to all zeroes, which produces a 3D ball. If you take a look at the mesh.obj that I put up (also see attached file below), you can make out part of the sphere, but the surface is all choppy. I'm not sure if I'm setting the dist variable or the normal variable correctly. Anyone have any suggestions?

Attached file:

mesh.obj

void addFractal(Field& field)
{
    float x_grid_max = 1.5;
    float y_grid_max = 1.5;
    float z_grid_max = 1.5;
    float x_grid_min = -x_grid_max;
    float y_grid_min = -y_grid_max;
    float z_grid_min = -z_grid_max;
    
    float z_w = 0;
    quaternion C;
    C.x = 0.0f;
    C.y = 0.0f;
    C.z = 0.0f;
    C.w = 0.0f;
    unsigned short int max_iterations = 8;
    float threshold = 4;
    
    string error_string;
    quaternion_julia_set_equation_parser eqparser;
    if (false == eqparser.setup("Z = sin(Z) + C * sin(Z)", error_string, C))
    {
        cout << "Equation error: " << error_string << endl;
        return;
    }
    else
    {
        cout << "Equation ok" << endl;
    }
    
    quaternion Z(x_grid_min, y_grid_min, z_grid_min, z_w);
    
	foreach3D(field.size(), [&](Vec3u upos)
    {
		auto& cell = field[upos];

        Z.x = upos.x / float(FieldSize - 1);
        Z.y = upos.y / float(FieldSize - 1);
        Z.z = upos.z / float(FieldSize - 1);
        
        cell.dist = eqparser.iterate(Z, max_iterations, threshold) - threshold;

        Z.x -= 0.01;
        float normal_x1 = eqparser.iterate(Z, max_iterations, threshold);
            
        Z.x += 0.02;
        float normal_x2 = eqparser.iterate(Z, max_iterations, threshold);
            
        Z.x -= 0.01;
        float normal_x = normal_x1 - normal_x2;
            
        Z.y -= 0.01;
        float normal_y1 = eqparser.iterate(Z, max_iterations, threshold);
            
        Z.y += 0.02;
        float normal_y2 = eqparser.iterate(Z, max_iterations, threshold);
            
        Z.y -= 0.01;
        float normal_y = normal_y1 - normal_y2;
            
        Z.z -= 0.01;
        float normal_z1 = eqparser.iterate(Z, max_iterations, threshold);
            
        Z.z += 0.02;
        float normal_z2 = eqparser.iterate(Z, max_iterations, threshold);
            
        Z.z -= 0.01;
        float normal_z = normal_z1 - normal_z2;
            
        normal_x /= 2.0;
        normal_y /= 2.0;
        normal_z /= 2.0;
        
        cell.normal = normalized(Vec3(normal_x, normal_y, normal_z));
    });
}