How Will I Compute The Normals at Vertices using 3D Object?

Pages: 12
Hello Professionals,

Good day. I would like to ask about the 3D aspect of my previous post recently. http://www.cplusplus.com/forum/general/233119/

Recently, I asked how to find the normals at vertices of a certain 2D object in order to set an offset not too far away from the boundary constraint and I am very honored to receive the solution from @lastchance.

Basically, I still have the same given values:

Given:
RED - boundary constraints
GREEN - normal constraints
BLUE - surface normal
d = distance between boundary constraints and normal constraints

I tried to create a sample picture because I cannot see a good example from internet.

Image: https://i.imgur.com/rvdZoRa.jpg

I would like to ask if I can still reuse the algorithm "behind" and "in front" algorithm in the 3D scenario, or do I have to change it? The previous code created by Mr. @lastchance was something like this:

Courtesy of Mr @lastchance:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
#include <iostream>
#include <iomanip>
#include <fstream>
#include <vector>
#include <cmath>
using namespace std;

//======================================================================

struct Pt2d { double x, y; };

Pt2d operator +( Pt2d p, Pt2d q )     { return { p.x + q.x, p.y + q.y }; }
Pt2d operator -( Pt2d p, Pt2d q )     { return { p.x - q.x, p.y - q.y }; }
Pt2d operator -( Pt2d p )             { return { -p.x     , -p.y      }; }
Pt2d operator *( double r, Pt2d p )   { return { r * p.x  , r * p.y   }; }
Pt2d operator /( Pt2d p, double r )   { return { p.x / r  , p.y / r   }; }
double dot( Pt2d p, Pt2d q ) { return p.x * q.x + p.y * q.y; }
double abs( Pt2d p ) { return sqrt( dot( p, p ) ); }

ostream &operator <<( ostream &strm, Pt2d p )
{
   const int w = 12;
   strm << fixed;
   return strm << setw( w ) << p.x << "  " << setw( w ) << p.y;
}

//======================================================================

vector<Pt2d> readData( string filename )
{
   vector<Pt2d> result;
   double x, y;
   ifstream fin( filename );
   while ( fin >> x >> y ) result.push_back( { x, y } );
   fin.close();
   return result;
}

//======================================================================

vector<Pt2d> setNormal( const vector<Pt2d> &boundary, bool anticlockwise )
{
   int N = boundary.size();
   vector<Pt2d> normal( N );

   for ( int i = 0; i < N; i++ )
   {
      int im = i - 1, ip = i + 1;                // Points behind and in front
      if ( i == 0     ) im = N - 1;              // Special cases at ends (assumes closed curve)
      if ( i == N - 1 ) ip = 0;

      // Get tangent (weighted on distances to nearby points)
      Pt2d t1 = boundary[ip]-boundary[i ];
      Pt2d t2 = boundary[i ]-boundary[im];
      double abst1 = abs( t1 );
      double abst2 = abs( t2 );
      Pt2d tangent = ( abst2 * t1 + abst1 * t2 ) / ( abst1 + abst2 ); 

      // Get normal vector by rotating 90 degrees
      Pt2d nvec = { tangent.y, -tangent.x };
      if ( !anticlockwise ) nvec = -nvec;
      normal[i] = nvec / abs( nvec );            // normalise
   }

   return normal;
}

//======================================================================

vector<Pt2d> setNormalConstraint( const vector<Pt2d> &boundary, const vector<Pt2d> &normal, double d )
{
   int N = boundary.size();
   vector<Pt2d> constraint( N );

   for ( int i = 0; i < N; i++ ) constraint[i] = boundary[i] - d * normal[i];
   return constraint;
}

//======================================================================

int main()
{
   bool anticlockwise = true;
   double d = 0.25;

   vector<Pt2d> boundary   = readData( "boundary.txt" );
   vector<Pt2d> normal     = setNormal( boundary, anticlockwise );
   vector<Pt2d> constraint = setNormalConstraint( boundary, normal, d );


   // Output
   #define SPACE << "      " <<
   ofstream fout( "output.txt" );
   for ( int i = 0; i < boundary.size(); i++ ) fout << boundary[i] SPACE normal[i] SPACE constraint[i] << '\n';
   fout.close();
}


Thank you for your valuable inputs in advance and God bless.


Last edited on
@kindgnice,
Can I ask if you really need the normals at the VERTICES, rather than the CENTROIDS of the triangles? What are you going to do with them? Any shading or texturing is likely to be done on the orientation of the flat triangle faces, whereas the vertices are pointy, and orientation isn't that well defined, as there may be an arbitrary number of triangles meeting at a vertex.

For a triangle you can find its vector area, and hence a normal vector, as half the vector cross product of two sides. The orientation depends on which way you ordered the points of the triangle. However, this is effectively located at the centroid of the triangle, not a vertex.

The only realistic way of getting normals at the vertices, rather than centroids, is some weighted average of the triangles meeting at that point.

The other thing you will have to worry about is the data structure giving points and triangles. In 2d a curve is just traced sequentially. There is no such concept of points "in front" or "behind" in 3d.

Your bunny figure is very complex. Why not start with a simpler polyhedron like tetrahedron, octahedron or icosahedron? They have a small number of triangular faces.
Thank you for your response @lastchance.

Can I ask if you really need the normals at the VERTICES, rather than the CENTROIDS of the triangles?


Oh yes, because this way, I can use the normals at vertices (i.e. vertices of 3D object), and make an offset not too far away from the coordinate (taking into consideration the normals at that particular coordinate/vertex.)

What are you going to do with them?


Taking the normals at vertex will help me generate the offset value not too far away from the coordinate/vertex with the help of normal. Generating these offset values will help me interpolate a surface between the offset value and the boundary constraint(vertex/coordinate). If I will take the centroid, it is not aligned(I am not sure if this is the proper term) to the coordinate/vertex.

The orientation depends on which way you ordered the points of the triangle.


Ohhh, how will I know the orientation/order? Basically I am taking OBJ samples freely online...

The only realistic way of getting normals at the vertices, rather than centroids, is some weighted average of the triangles meeting at that point.


Ah, for 2D, weighted average distances at a certain point and for 3D, weighted average areas at a certain point? Hmm, this sounds complicated for me, because upon observing different OBJ files before, some points contain few triangles, some points contain much more triangles coinciding to a specific point...

Your bunny figure is very complex. Why not start with a simpler polyhedron like tetrahedron, octahedron or icosahedron? They have a small number of triangular faces.


Oh yes, I agree with you. Probably, an icosahedron would be a good start...
https://github.com/nopjia/tracer/blob/master/data/icosahedron.obj

@kindgnice,
Please give a clearer picture of what you are actually trying to do AT THE END OF THE DAY.

Are you simply trying to put a fixed-thickness border inside your object? If so, then this has been done in this forum before in 2d:
http://www.cplusplus.com/forum/beginner/223356/#msg1023629
For that you DO NOT NEED TO FIND NORMALS AT ALL.

That approach does not generalise to 3d if more than three triangles meet at a point, however. So, it is vitally important that you tell us what your END OBJECTIVE is, NOT how you THINK you should approach it.


FWIW, below is the border polygon code for your (2d) bunny. The output is
https://imgur.com/a/MvuoC1g

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
#include <iostream>
#include <iomanip>
#include <fstream>
#include <vector>
#include <cmath>
#include <string>
#include <cstdlib>
using namespace std;

const double SMALL = 1.0e-10;

// Structures/classes
struct Point{ double x, y; };

// Prototypes
vector<Point> readBoundary( istream &in );
vector<Point> border( const vector<Point> &p, double thickness );
double signedArea( const vector<Point> &p );
void output( ostream &strm, const vector<Point> &p, const vector<Point> &q );


//======================================================================


int main()
{
   ifstream in( "input.txt" );                   // Put your polygon coordinates in here
   double thickness = 0.02;                      // Border thickness - ADJUST FOR YOUR PARTICULAR PROBLEM

   vector<Point> p = readBoundary( in );         // Read vertices
   vector<Point> q = border( p, thickness );     // Compute internal vertices

   ofstream out( "output.txt" );                 // Output to file
   output(  out, p, q );
   out.close();

   cout << "\nNumber of points: " << p.size() << '\n';
   cout << "\nEnclosed area (-ve means points are clockwise): " << signedArea( p ) << '\n';
}


//======================================================================


vector<Point> readBoundary( istream &in )
{
   vector<Point> p;
   Point pt;

   while ( in >> pt.x >> pt.y ) p.push_back( pt );         // Read points as x y
   if ( p.size() == 0 )                                    // Deal with unable-to-read-points case
   {
      cout << "No points";
      exit( 1 );
   }

   pt = p.back();                                          
   if ( pt.x != p[0].x || pt.y != p[0].y ) p.push_back( p[0] );      // Make sure that polygon is closed

   return p;
}


//====================================================================== 


vector<Point> border( const vector<Point> &p, double thickness )
{
   //=====================================================//
   // Assumes:                                            //
   //    N+1 vertices:   p[0], p[1], ... , p[N-1], p[N]   //
   //    Closed polygon: p[0] = p[N]                      //
   //    No zero-length sides                             //
   // Returns:                                            //
   //    Internal poly:  q[0], q[1], ... , q[N-1], q[N]   //
   //=====================================================//
   int N = p.size() - 1;
   vector<Point> q(N+1);              

   double a, b, A, B, d, cross;

   double displacement = thickness;
   if ( signedArea( p ) < 0 ) displacement = -displacement;     // Detects clockwise order

   // Unit vector (a,b) along last edge
   a = p[N].x - p[N-1].x;
   b = p[N].y - p[N-1].y;
   d = sqrt( a * a + b * b );  
   a /= d;
   b /= d; 

   // Loop round the polygon, dealing with successive intersections of lines
   for ( int i = 0; i < N; i++ )
   {
      // Unit vector (A,B) along previous edge
      A = a;
      B = b;
      // Unit vector (a,b) along next edge
      a = p[i+1].x - p[i].x;
      b = p[i+1].y - p[i].y;
      d = sqrt( a * a + b * b );
      a /= d;
      b /= d;
      // New vertex
      cross = A * b - a * B;
      if ( abs( cross ) < SMALL )      // Degenerate cases: 0 or 180 degrees at vertex
      {
         q[i].x = p[i].x - displacement * b;
         q[i].y = p[i].y + displacement * a;
      }
      else                             // Usual case
      {
         q[i].x = p[i].x + displacement * ( a - A ) / cross;
         q[i].y = p[i].y + displacement * ( b - B ) / cross;
      }
   }

   // Close the inside polygon
   q[N] = q[0];

   return q;
}


//======================================================================


double signedArea( const vector<Point> &p )
{
   //========================================================//
   // Assumes:                                               //
   //    N+1 vertices:   p[0], p[1], ... , p[N-1], p[N]      //
   //    Closed polygon: p[0] = p[N]                         //
   // Returns:                                               //
   //    Signed area: +ve if anticlockwise, -ve if clockwise //
   //========================================================//
   int N = p.size() - 1;
   double A = 0;
   for ( int i = 0; i < N; i++ ) A += p[i].x * p[i+1].y - p[i+1].x * p[i].y;
   A *= 0.5;
   return A;
}


//======================================================================


void output( ostream &strm, const vector<Point> &p, const vector<Point> &q )
{
   int N = p.size() - 1;

   // Formatting (adapt to your needs)
   strm.setf( ios::fixed );                 // Fixed or scientific
   strm.precision( 6 );                     // Decimal digits after point
   #define SP << " " << setw( 12 ) <<       // Spacer and field width

   for ( int i = 0; i <= N; i++ ) strm SP p[i].x SP p[i].y SP q[i].x SP q[i].y << '\n';
}


//======================================================================  



input.txt (minus point which is clearly wrong)
0.2137          0.386
0.177           0.386
0.1437          0.4016
0.11            0.418
0.078           0.443
0.0485          0.474
0.047           0.519
0.061           0.562
0.085           0.604
0.108           0.643
0.1345          0.6794
0.165           0.7126
0.201           0.73
0.235           0.749
0.232           0.798
0.24            0.839
0.252           0.8775
0.27            0.921
0.2984          0.9575
0.3303          0.9815
0.369           0.984
0.3896          0.9466
0.381           0.9042
0.3655          0.8713
0.3506          0.832
0.3464          0.7894
0.375           0.8243
0.392           0.8586
0.417           0.8915
0.4497          0.9142
0.4896          0.9124
0.512           0.8786
0.503           0.833
0.4812          0.789
0.4504          0.758
0.4155          0.7349
0.3808          0.713
0.35            0.6884
0.3697          0.6576
0.39            0.6285
0.42            0.6036
0.4514          0.63
0.486           0.649
0.521           0.659
0.5544          0.669
0.589           0.6682
0.6213          0.662
0.655           0.654
0.689           0.6425
0.7196          0.6238
0.7545          0.6046
0.7853          0.576
0.8182          0.5407
0.84            0.5068
0.8605          0.4698
0.8717          0.4243
0.879           0.3825
0.8762          0.3385
0.869           0.2903
0.8995          0.3095
0.9332          0.3032
0.955           0.2675
0.9617          0.229
0.9526          0.1816
0.926           0.146
0.8974          0.1194
0.8647          0.117
0.837           0.138
0.7998          0.1233
0.7886          0.0866
0.762           0.0533
0.7286          0.0472
0.6946          0.0448
0.6612          0.041
0.6292          0.0384
0.5935          0.038
0.564           0.0283
0.5313          0.039
0.5327          0.0843
0.545           0.123
0.513           0.1227
0.4804          0.125
0.4467          0.14
0.4106          0.148
0.3868          0.1157
0.379           0.0767
0.3564          0.0373
0.3213          0.032
0.2883          0.0367
0.2514          0.049
0.2433          0.0915
0.2706          0.128
0.3045          0.145
0.3026          0.189
0.27            0.2114
0.251           0.253
0.241           0.297
0.23            0.342

Last edited on
OK @kindgnice, I had another look at this. This is a 2-part post because of the length.

I don't use OpenGL (yet) but I've been skimming through the .obj format at
https://en.wikipedia.org/wiki/Wavefront_.obj_file

It seems there is scope there for outputting "vertex normals". These are the lines starting "vn". Indeed, the sample file of an icosahedron that you pointed me to allegedly includes them. The only problem is that they are clearly NOT VERTEX normals, as an icosahedron has 20 faces and 12 vertices: in fact, they are the 20 FACE normals (as computed by my code in the following post).

So, in the post that follows I include code to read an .obj file (called input.obj). This is the file for an icosahedron. It then computes both face area normals (which I have checked against the ones in the file) and vertex normals (which I have computed by averaging the normals for triangular faces meeting at each vertex, although there is an option to weight them more appropriately on enclosed angle if desired). At present this is written to cout (the console) - you can redirect the ostream to a file if you want to. See the code (near the end of main() ).

input.obj
# icosahedron.obj

g icosahedron

v 0.000000 -0.525731 0.850651
v 0.850651 0.000000 0.525731
v 0.850651 0.000000 -0.525731
v -0.850651 0.000000 -0.525731
v -0.850651 0.000000 0.525731
v -0.525731 0.850651 0.000000
v 0.525731 0.850651 0.000000
v 0.525731 -0.850651 0.000000
v -0.525731 -0.850651 0.000000
v 0.000000 -0.525731 -0.850651
v 0.000000 0.525731 -0.850651
v 0.000000 0.525731 0.850651

vn 0.934172 0.356822 0.000000
vn 0.934172 -0.356822 0.000000
vn -0.934172 0.356822 0.000000
vn -0.934172 -0.356822 0.000000
vn 0.000000 0.934172 0.356822
vn 0.000000 0.934172 -0.356822
vn 0.356822 0.000000 -0.934172
vn -0.356822 0.000000 -0.934172
vn 0.000000 -0.934172 -0.356822
vn 0.000000 -0.934172 0.356822
vn 0.356822 0.000000 0.934172
vn -0.356822 0.000000 0.934172
vn 0.577350 0.577350 -0.577350
vn 0.577350 0.577350 0.577350
vn -0.577350 0.577350 -0.577350
vn -0.577350 0.577350 0.577350
vn 0.577350 -0.577350 -0.577350
vn 0.577350 -0.577350 0.577350
vn -0.577350 -0.577350 -0.577350
vn -0.577350 -0.577350 0.577350

f 2//1 3//1 7//1
f 2//2 8//2 3//2
f 4//3 5//3 6//3
f 5//4 4//4 9//4
f 7//5 6//5 12//5
f 6//6 7//6 11//6
f 10//7 11//7 3//7
f 11//8 10//8 4//8
f 8//9 9//9 10//9
f 9//10 8//10 1//10
f 12//11 1//11 2//11
f 1//12 12//12 5//12
f 7//13 3//13 11//13
f 2//14 7//14 12//14
f 4//15 6//15 11//15
f 6//16 5//16 12//16
f 3//17 8//17 10//17
f 8//18 2//18 1//18
f 4//19 10//19 9//19
f 5//20 9//20 1//20





The output is the following. In the final section, the six floating-point columns are the vertices components and the vertex-unit-normal components. In this particular example ... they are the same!! This is because of the symmetry (and radius 1) of the given icosahedron.
Note that:
(1) the vertex indices for triangles have been reduced by 1 (because of counting array elements from 0).
(2) the normals have unit length (i.e. 1.0). If you want to plot them for an object with much smaller dimensions then you will have to scale them before plotting.

Vertices: 12
0: 0  -0.525731  0.850651
1: 0.850651  0  0.525731
2: 0.850651  0  -0.525731
3: -0.850651  0  -0.525731
4: -0.850651  0  0.525731
5: -0.525731  0.850651  0
6: 0.525731  0.850651  0
7: 0.525731  -0.850651  0
8: -0.525731  -0.850651  0
9: 0  -0.525731  -0.850651
10: 0  0.525731  -0.850651
11: 0  0.525731  0.850651

Triangles: 20
0: 1 2 6 0.934172  0.356822  0
1: 1 7 2 0.934172  -0.356822  0
2: 3 4 5 -0.934172  0.356822  0
3: 4 3 8 -0.934172  -0.356822  -0
4: 6 5 11 0  0.934172  0.356822
5: 5 6 10 0  0.934172  -0.356822
6: 9 10 2 0.356822  0  -0.934172
7: 10 9 3 -0.356822  -0  -0.934172
8: 7 8 9 -0  -0.934172  -0.356822
9: 8 7 0 0  -0.934172  0.356822
10: 11 0 1 0.356822  0  0.934172
11: 0 11 4 -0.356822  0  0.934172
12: 6 2 10 0.57735  0.57735  -0.57735
13: 1 6 11 0.57735  0.57735  0.57735
14: 3 5 10 -0.57735  0.57735  -0.57735
15: 5 4 11 -0.57735  0.57735  0.57735
16: 2 7 9 0.57735  -0.57735  -0.57735
17: 7 1 0 0.57735  -0.57735  0.57735
18: 3 9 8 -0.57735  -0.57735  -0.57735
19: 4 8 0 -0.57735  -0.57735  0.57735

Vertices and vertex normals:
            0     0.000000    -0.525731     0.850651     0.000000    -0.525731     0.850651
            1     0.850651     0.000000     0.525731     0.850651     0.000000     0.525731
            2     0.850651     0.000000    -0.525731     0.850651     0.000000    -0.525731
            3    -0.850651     0.000000    -0.525731    -0.850651     0.000000    -0.525731
            4    -0.850651     0.000000     0.525731    -0.850651     0.000000     0.525731
            5    -0.525731     0.850651     0.000000    -0.525731     0.850651     0.000000
            6     0.525731     0.850651     0.000000     0.525731     0.850651     0.000000
            7     0.525731    -0.850651     0.000000     0.525731    -0.850651     0.000000
            8    -0.525731    -0.850651     0.000000    -0.525731    -0.850651     0.000000
            9     0.000000    -0.525731    -0.850651     0.000000    -0.525731    -0.850651
           10     0.000000     0.525731    -0.850651     0.000000     0.525731    -0.850651
           11     0.000000     0.525731     0.850651     0.000000     0.525731     0.850651

Last edited on
SECOND OF TWO-PART POST - CODE FOR THE ABOVE
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
#include <iostream>
#include <fstream>
#include <sstream>
#include <iomanip>
#include <vector>
#include <cmath>
#include <string>
using namespace std;

const double PI = 3.141592653589793;
const bool ANTICLOCKWISE = true;                           // Traversal of triangles, seen from OUTSIDE

//================================================

struct Pt                                                  // Coordinates of a vector or point
{
    double x, y, z;
};                              

struct Triangle                                            // Holds the INDICES and vector area of triangles
{ 
   int v[3];
   Pt area;
   Pt normal;
};

Pt operator +( Pt p, Pt q ) { return { p.x + q.x, p.y + q.y, p.z + q.z }; }
Pt operator -( Pt p, Pt q ) { return { p.x - q.x, p.y - q.y, p.z - q.z }; }
Pt operator *( double r, Pt p ) { return { r * p.x, r * p.y, r * p.z }; }
Pt operator *( Pt p, double r ) { return r * p;                         }
Pt operator /( Pt p, double r ) { return { p.x / r, p.y / r, p.z / r }; }
double dot( Pt p, Pt q ) { return p.x * q.x + p.y * q.y + p.z * q.z; }
double abs( Pt p ) { return sqrt( dot( p, p ) ); }
Pt cross( Pt p, Pt q ) { return { p.y * q.z - p.z * q.y, p.z * q.x - p.x * q.z, p.x * q.y - p.y * q.x }; }

ostream &operator << ( ostream &strm, Pt  p ) { return strm << p.x << "  " << p.y << "  " << p.z; }
istream &operator >> ( istream &strm, Pt &p ) { return strm >> p.x >> p.y >> p.z; }

//================================================


void readFile( string filename, vector<Pt> &vertices, vector<Triangle> &triangles )
{
   vertices.clear();                                       // Comment out if adding to an existing collection
   triangles.clear();                                      //

   string line;

   ifstream in( filename );
   while ( getline( in, line ) )                           // Read whole line
   {
      if ( line.find_first_of( "vVfF" ) == string::npos ) continue;     // skip pointless lines

      istringstream ss( line );                            // Put line into a stream for input
      string s;
      ss >> s;                                             // Get identifier
      if ( s == "v" || s == "V" )                          // VERTICES
      {
         Pt p;
         ss >> p;
         vertices.push_back( p );
      }
      else if ( s == "f" || s == "F" )                     // FACES
      {
         string si, sj, sk;
         ss >> si >> sj >> sk;                             // FIRST, read into individual strings
         Triangle t;
         t.v[0] = stoi( si ) - 1;                          // Get the FIRST integer from each string
         t.v[1] = stoi( sj ) - 1;                          // NOTE: subtract 1 because arrays start from 0
         t.v[2] = stoi( sk ) - 1;
         triangles.push_back( t );
      }
   }
   in.close();
}


//======================================================================


void calcFaces( const vector<Pt> &vertices, vector<Triangle> &triangles )
{
   for ( int t = 0; t < triangles.size(); t++ )
   {
      int v0 = triangles[t].v[0];                          // Global vertex index
      int v1 = triangles[t].v[1];
      int v2 = triangles[t].v[2];
      Pt side1 = vertices[v1] - vertices[v0];              // Side vectors meeting at this vertex
      Pt side2 = vertices[v2] - vertices[v0];
      triangles[t].area = 0.5 * cross( side1, side2 );     // Triangle vector area is half the cross product of the sides
      if ( !ANTICLOCKWISE ) triangles[t].area = -1.0 * triangles[t].area;
      double A = abs( triangles[t].area );
      triangles[t].normal = triangles[t].area / A;         // Unit normal vector
   }
}


//================================================


vector<Pt> calcNormals( const vector<Pt> &vertices, const vector<Triangle> &triangles )
{
   vector<Pt> vnormal( vertices.size(), { 0.0, 0.0, 0.0 } );

   for ( int t = 0; t < triangles.size(); t++ )            // Loop through triangle list
   {
      Pt norm = triangles[t].area;
      norm = norm / abs( norm );                           // Unit normal vector for this triangle

      for ( int i = 0; i < 3; i++ )                        // Update weighted sum at each vertex
      {
         int v0 = triangles[t].v[i];                       // Global vertex index
         int v1 = triangles[t].v[(i+1)%3];                 
         int v2 = triangles[t].v[(i+2)%3];
         Pt side1 = vertices[v1] - vertices[v0];           // Side vectors meeting at this vertex
         Pt side2 = vertices[v2] - vertices[v0];
//       double weight = acos( dot( side1, side2 ) / ( abs( side1 ) * abs( side2 ) ) );   // Angle enclosed at this vertex
         double weight = 1.0;                              // Use this if you want unweighted averages
         vnormal[v0] = vnormal[v0] + weight * norm;        // Contribution to normal vector
      }
   }

   for ( int v = 0; v < vertices.size(); v++ )
   {
      vnormal[v] = vnormal[v] / abs( vnormal[v] );         // Normalise to unit length
   }

   return vnormal;
}


//================================================


void output( ostream &out, const vector<Pt> &vertices, const vector<Pt> &vnormal )
{
   int N = vertices.size();

   // Formatting (adapt to your needs)
   out.setf( ios::fixed );                  // Fixed or scientific
   out.precision( 6 );                      // Decimal digits after point
   #define SP << " " << setw( 12 ) <<       // Spacer and field width

   for ( int i = 0; i < N; i++ )
   {
      out SP i
          SP vertices[i].x SP vertices[i].y SP vertices[i].z
          SP vnormal [i].x SP vnormal [i].y SP vnormal [i].z << '\n';
   }
}


//================================================


void summary( const vector<Pt> &vertices, const vector<Triangle> &triangles )
{
   cout << "\nVertices: " << vertices.size() << '\n';
   for ( int i = 0; i < vertices.size(); i++ ) 
   {
      cout << i << ": " << vertices[i] << '\n';
   }

   cout << "\nTriangles: " << triangles.size() << '\n';
   for ( int t = 0; t < triangles.size(); t++ ) 
   {
      cout << t << ": ";
      for ( int i = 0; i < 3; i++ )
      {
         cout << triangles[t].v[i] << " ";
      }
      cout << triangles[t].normal << '\n';
   }
}


//================================================


int main()
{
   string filename = "input.obj";
   vector<Pt> vertices;
   vector<Triangle> triangles;

   readFile( filename, vertices, triangles );                   // Read vertices and triangles
   calcFaces( vertices, triangles );                            // Calculate area vectors for triangles
   summary( vertices, triangles );                              // Check

   vector<Pt> vnormal = calcNormals( vertices, triangles );     // Calculate normals at vertices
   cout << "\nVertices and vertex normals:\n";
   output( cout, vertices, vnormal );
}

//================================================  

Last edited on
icosahedron that you pointed me to allegedly includes them. The only problem is that they are clearly NOT VERTEX normals, as an icosahedron has 20 faces and 12 vertices: in fact, they are the 20 FACE normals

Let me add to that.

That icosahedron has clearly its center of mass at (0,0,0) to keep things simple. Furthermore, its radius is 1.
A vertex is essentially a 3D vector from (0,0,0) to the position of the vertex. With radius 1 the length of those vectors is 1.

A normal is a 3D vector too, but its length should be kept 1. A normal is a direction. It has nothing to do with position.

With radius 1 the vertices and vertex normals look identical.

If we move the icosahedron "left" by two units, we add to each vertex vector (-2,0,0). The normals will not change.

If we want to increase radius to 3, we would update the vertices of the moved object with
v' = 3 * ( v - {-2,0,0} ) + {-2,0,0}


When you render (draw) a pixel of a triangle, the color of that pixel depends on how light hits that pixel. Surface reflects light. Reflection depends on the angle of impact. We can compute the angle from the normal of the surface and the direction of light. In real life the normal is orthogonal to the surface.

A face of icosahedron is a triangle. Triangle is planar surface. Flat. In principle every pixel of it has same normal. The FACE normal.

You can give OpenGL (or other rendering system) three (vertex,normal) pairs and then it can draw a triangle. You can give same face normal in each pair. When the system considers a pixel on the triangle, it will interpolate the position from the three vertices and normal from the three normals. If all three normals are identical, then every pixel has the same normal, the face normal.

It is quite easy to see each triangle from the image because two triangles that share edge can have quite different normals.

If you want "boundary" at 0.2 "above" triangle, you simply add 0.2*face_normal to the vertices of the triangle. You simply move triangle "up".
If you do the same for the neighboring triangle, it will move too, but there will be a gap between the edges because the triangles move to different directions.


A vertex can be member of multiple triangles. It will have a different normal, depending on which face we are looking at.


Now the vertex normal. There is only one per vertex. All vertices of a triangle can have different normal. All pixels of the triangle will have different normals and different color. Edges between triangles become less clear. The icosahedron will look more smooth, more round, even though it still is made of "flat" triangles.

If you move each vertex by 0.2*vertex_normal, the triangle will probably move less that 0.2, but it will stretch out to larger size and will continue to share vertices and edges with its neighbours.
Thank you for this, @Keskiverto. I think your two paragraphs sum up the need for vertex normals as well as face normals. (The emphasised bits below are mine.)

(1) If you want "boundary" at 0.2 "above" triangle, you simply add 0.2*face_normal to the vertices of the triangle. You simply move triangle "up".
If you do the same for the neighboring triangle, it will move too, but there will be a gap between the edges because the triangles move to different directions.


(2) If you move each vertex by 0.2*vertex_normal, the triangle will probably move less that 0.2, but it will stretch out to larger size and will continue to share vertices and edges with its neighbours.


I think if you only have three (non-coplanar) triangles meeting at a point then the three displaced and extended, but parallel, planes will continue to meet at a point ... however, if there are more than three triangles then they won't. As you say, there will be a gap and the nicely triangulated surface will fall to pieces. So you need vertex normals in 3d.

My only problem was how to define a "vertex normal" at what is, essentially, a local peak (or trough). In my code I simply added, vectorially, the normals from adjacent triangles, and then normalised to unit length. I commented out an alternative version with weightings based on enclosed angle, which might be better for badly-triangulated surfaces.

The face normals are also available from the triangles in the code. There is no ambiguity in those because the triangle faces are planar. Also available are the area vectors; (in the direction of the outward normal, but magnitude of the triangle area). These have a lot of uses, including finding the total volume of the figure as
(1/3)*sum(faces)A.r
where A is the area vector of the face and r any position vector on it (e.g. a vertex).
Last edited on
Thank you for your program codes and questions @lastchance. It took me some time to think how to explain my idea. Hmm, I found some of my 3D samples and I am not sure if this would give something intuitively about my situation...

I did try to "manually" set the 3D coordinate points by hand using his 3D heart shape shown here in this picture.

Original heart OBJ: https://i.imgur.com/x36RnDU.png

Since I don't know how to generate the vertex at normals, what I did was I tried to "scale down" (i.e. 0.99 scale or slightly "smaller" than the original OBJ shape) the heart OBJ which will give me approximately an offset distance from the original coordinates of the OBJ. Doing this so, it helped me to run interpolate the implicit surface something like this...

Output after scaling down(0.99 of the original): https://i.imgur.com/sH5wHxw.png

However, if at some point, I set the offset points (i.e. anywhere or not close to the coordinate points), my output shows like a destroyed version of heart: https://i.imgur.com/CDBcbH4.png


Last edited on
Well, my code should give you the normals at the vertices. They are written in routine output and you can adapt that to your needs. If you want displaced points then use (vectorially)
vertex - d* normal
where d is the distance that you want to move inward. Given the dimensions of your object, d should be small.
Thank you very much @keskiverto and @lastchance for your valuable inputs. I think I need to re-read and review the codes once again to digest and understand what's happening... :( But I am very honored and thankful and grateful to receive your guidance as always... I hope I can understand these codes soon...
I am so SORRY to everyone Sir @Keskiverto and Sir @Lastchance for my late response. It took me so long to respond due to several trials of the codes, presentations and other subjects. :'(
Apologies for that...

@lastchance: Did you mean Sir that if my 3D object is not symmetric, then the resulting vertex normals will be different from the vertex coordinates?

@keskiverto: Did you mean Sir that using this formula "0.2*vertex_normal" will help me scale the object?

Last edited on
Yes, you are right Sir @lastchance. After trying to run the program with heart 3D shape, the results are different (meaning the vertices values are different from the vertex normals).
Last edited on
Hello Sir @lastchance. I tried to incorporate your codes in my program and I would like to share my output. I got some outputs like these (using distance offset value as 0.5 units):

Original Heart OBJ: https://i.imgur.com/EK7uySj.png

Generated Heart OBJ (adding offsets formula [vertex-d*normal] from the code): https://i.imgur.com/BjX4XSe.png

Original Star OBJ: https://i.imgur.com/pHzAkbQ.png

Generated Star OBJ (adding offsets formula [vertex-d*normal] from the code): https://i.imgur.com/WwryFe8.png

Original Deer OBJ: https://i.imgur.com/uNudWYa.png

Generated Deer OBJ (adding offsets formula [vertex-d*normal] from the code):

None. (0 KB) :'(

1. I was wondering why it didn't work for other OBJs? Is there a specific reason why?
2. Is it because of the construction itself of a certain OBJ?
3. Does it have something to do with the normals? Does it sometimes find difficulty to decide whether the certain "offseted" constraint is inside or outside the surface boundary of OBJ?

Hmm, I mean, if I will use the formula vertex-d*normal, will I have to make an additional IF statement for the offset to check if the point is inside or outside the surface? Or the code already knows that the offseted point is inside the surface?

(Apologies I couldn't find a good example drawing from internet that's why I tried to draw my sample by myself using Sketchpad). For example if we are talking about complicated surface. Here's my sample picture:

https://i.imgur.com/YFC75jc.png

Does the code already know that the offset value I will compute using the formula vertex-d*normal is "inside" the surface (not outside), or I will have to add another condition Sir?

Thank you very much for your response in advance and God bless. :)


Last edited on
@kindgnice,
I can't debug your programs for you. If it is generating NO output for a specific case then something is wrong with your file handling. If you provide the input file I'll look at it.

The size of d depends on how big the original figure is. All the normals are UNIT normals: that is, they have length 1 unit. Choose a size for d to allow for the size of your figure, not the same for everything.

The direction of the normals depends on which order the triangles are traversed in. I'm following the convention specified for .OBJ files ... which doesn't necessarily mean that whoever compiled those files did. You can do a global change of direction simply by setting ANTICLOCKWISE to false.
Thank you for your response Sir @lastchance.

I tried to give different values of offset to different OBJ files.

(a) HEART
For the HEART obj, if I give 0.5 offset/distance value, here is the output or generated surface I got: https://i.imgur.com/GAryLPE.png

If I will give 0.1 offset/distance value, the output that I got has some unexplained structure at the bottom which I can't explain why :'( https://i.imgur.com/6KfhUN7.png

If I will give 0.01 offset/distance value, then the output I got will be almost a look alike of 0.5 offset/distance value generated: https://i.imgur.com/XAviDO6.png

(b) STAR
On the other hand, if I will use a different OBJ file, for example, a STAR obj, if I will give a 0.5 offset value, the output will be: https://i.imgur.com/RzP7SPR.png

If I will give 0.1 offset value, the output will be (I don't know Sir why it looks creepy): https://i.imgur.com/eLVKm07.png

And if I will use a 0.01 offset value, the output is even worst: https://i.imgur.com/TafQe2f.png

I am wondering why the STAR obj has an unexpected behavior when I reduce the offset value from 0.5 to 0.01. Honestly, according to the research paper I read, it says that 0.01 unit offset/distance between boundary constraint and normal constraint is a good choice for any OBJ files. Unfortunately, it didn't work well with STAR obj, and I was wondering why....
------------------------------------------------------------------------------------------------

*ANOTHER SCENARIO*

If I will use a DEER obj, if I try to give 0.5 offset value, there is no generated surface output. (e.g. 0 KB). Here's the OBJ file of the deer.

Picture: https://i.imgur.com/fydBpPk.png
OBJ file: https://drive.google.com/file/d/1NTQlG10aQFVkgBD-cfaL9yt3KKyt8dHi/view?usp=sharing

------------------------------------------------------------------------------------------------

*ANOTHER OBSERVATION*

I also observed Sir regarding the text details of OBJ files. (using distance/offset=0.5)
For example, if we are going to talk about the STAR obj, if I will try to open the .txt format of the original STAR obj file, the 3D output and details I see are as follows.

Output: https://i.imgur.com/HLtTOYf.png [I don't know why it behaves like that Sir?]
OBJ file (filename: star.obj): https://drive.google.com/file/d/1Y_ZeHbt2a-43UMbCLuXcDuWI4aSQFLF8/view?usp=sharing

On the other hand, when I tried to remove the "vt (textures)" and "vn (normals)" from the original STAR obj file using MeshLab, what happens is that, the text details remained in the obj file were only "v" and "f". In this case, when I tried to run this new extracted OBJ file, the output I got is (which is the same OBJ I used and showed from the figure above for 0.5 offset/distance for letter (b) scenario):

Output: https://i.imgur.com/RzP7SPR.png
OBJ file (filename: star_using_MeshLab.obj): https://drive.google.com/file/d/1hfiZcTjmr5RDlRUMP_EvKZAQJ43LY3xD/view?usp=sharing

------------------------------------------------------------------------------------------------
I also tried to change the ANTICLOCKWISE state to false Sir, but I think the output is the same as I had before (I add "-" sign to the distance if I changed the ANTICLOCKWISE to false). I also did try to delete the "vt" and "vn" from the DEER obj using MeshLab to see if there's a hope to interpolate surface like what happened to star_using_MeshLab.obj but unfortunately, there was no output as well (0 KB).



Last edited on
@kindgnice,
Are you sure that you are using the routines from the latest post here:
http://www.cplusplus.com/forum/general/235547/#msg1056101
and NOT mixed in earlier versions of various routines.

In particular, be careful with the readFile routine. In the latest version this read the first item into a string and hence should easily be able to distinguish between "v", "vt" and "vn":
1
2
3
4
5
6
7
8
9
10
11
12
13
   while ( getline( in, line ) )                           // Read whole line
   {
      if ( line.find_first_of( "vVfF" ) == string::npos ) continue;     // skip pointless lines

      istringstream ss( line );                            // Put line into a stream for input
      string s;     //<=====
      ss >> s;                                             // Get identifier
      if ( s == "v" || s == "V" )                          // VERTICES
      {
         Pt p;
         ss >> p;
         vertices.push_back( p );
      }


Earlier versions (where I was not aware of the alternate vt and vn items) just picked out the first element as a char (single character). That would explain why some files produced zero output and why some files worked only when vt and vn items were removed.

DO NOT MIX AND MATCH DIFFERENT VERSIONS - THEY ARE DESIGNED TO WORK AS A WHOLE.



If a file contains "vn" items then these "should" be the vertex normals anyway: you shouldn't need my routines to calculate them.


I repeat my earlier point that a single value of d is NOT going to work for all files - it depends on how big your figures are.

Also, your whole strategy of going inward or outward a fixed distance is NOT a good way of SCALING a figure, as it is BOUND to distort its shape! IF YOU SIMPLY WANT TO SCALE A FIGURE THEN SIMPLY MULTIPLY ITS VERTEX VECTORS BY THE SCALE FACTOR! Things like normal vectors are for playing with lighting, viewpoints etc, NOT for scaling.
Last edited on
Are you sure that you are using the routines from the latest post here:
http://www.cplusplus.com/forum/general/235547/#msg1056101
and NOT mixed in earlier versions of various routines.


Yes Sir, I rechecked the readFile routine and the ones that I used was as described in the link you recently provided and not the ones from the 2D-related program before.

Also, your whole strategy of going inward or outward a fixed distance is NOT a good way of SCALING a figure, as it is BOUND to distort its shape!


Honestly Sir, I had asked this question regarding the inside or outside (https://i.imgur.com/YFC75jc.png) offseted normal constraint because I am not sure Sir if the routine outputs the inside or outside constraint by using the formula vertex-d*normal.

I would also like to share Sir one of my observations using another OBJ file. It is an CAT obj file. The obj file is: https://drive.google.com/file/d/13B518o27HFK84pWIapQZvckFsXM8ZqT-/view?usp=sharing
What happens Sir is that, when I try to run the program, in the normals, there is/are "-nan(ind)" outputs.. Such as shown below:


The text details are too long Sir. I will try to share the output via this link:
https://justpaste.it/250uj


And the computed offset normal constraints had the "-nan(ind)" values too.


The text details are too long Sir. I will try to share the output via this link:
https://justpaste.it/1yvi1


I am not sure Sir if these "-nan(ind)" values are one of the reasons why I couldn't get outputs... :(

I am also speculating Sir if there are problems if the OBJ files have "sharp" edges. Like for example in this figure: https://i.imgur.com/0ff6XYA.png

I am speculating maybe in some cases if the surface of an object contains very sharp edges like this, maybe the interpolation fails? Because as the computed normals are getting closer with the other normals, the charge or offset values tend to “repel” at some point which makes other portion of the figure explode?

For example, if I set the offset value a small value/charge such as 0.1 [inside where f(x)>0], for instance, it will come to a point where if the edges like as shown in the left figure comes closer, the charges will “repel” and get confused? Unlike the other OBJ files I used which has no sharp edges like HEART? I am also not sure Sir if my speculation has meaning… :(

Note: The inside surface is f(x)>0 and outside surface is f(x)<0. Surface is f(x)=0
Last edited on
kindgnice wrote:
I would also like to share Sir one of my observations using another OBJ file. It is an CAT obj file. The obj file is: https://drive.google.com/file/d/13B518o27HFK84pWIapQZvckFsXM8ZqT-/view?usp=sharing
What happens Sir is that, when I try to run the program, in the normals, there is/are "-nan(ind)" outputs..


Why don't you look at this particular input file. Then you will see why. What shape are the faces if they have four vertices?

The code as written assumes there are just three vertices following an "f" for face; i.e. faces are triangles. However, most of the faces in that particular input file have four vertices; (actually, some have three and some have five). Hence, some of the vertices are never assigned to any triangles and consequently the code as currently written will not find normals at them.

Sorry, but I'm not rewriting my code to deal with arbitrary polygons as faces. Over to you. If you want to use these type of cases then you can break the polygonal faces into triangles.

At a rough guess, this is probably happening with some of your other cases.


Your deer file is spurious. It has, for example, two triangles with vertex triplets 480, 249, 74 and 480, 74, 249. This produces two vector areas pointing in opposite directions, which add up to give a vector normal of zero. Maybe this is some 2-d part of your picture, but the code isn't designed to work with this.


I'm sorry, but it's up to you to fix your input files. Make sure that all faces are triangles and the figures are genuinely three-dimensional.
Last edited on
Thank you very much Sir for your always support. I have been trying to look for methods on how to deal with the many faces of each OBJ file(s), and I have found the software which can do the trick of making all the faces stick to 3 faces only. And that was Blender.

However, when I succeeded to get the output 3 faces for the same OBJ file, and tried running the codes again, there was no interpolation surface. You are also right regarding the selection of proper distance/offset and you are also correct that this value differs for each OBJ file. My problem is how do I select the proper distance/offset. There were times I was just guessing (like trial and error). Sometimes I get the right offset, most of the time, I have to give different values like many attempts of trials and errors before I arrived with the right offset value...
Pages: 12