diff --git a/entity.cpp b/entity.cpp index b55ea84..7933439 100644 --- a/entity.cpp +++ b/entity.cpp @@ -40,22 +40,9 @@ bool Entity::HasPlane(void) { void Entity::PlaneGetExprs(ExprVector *n, Expr **dn) { if(type == WORKPLANE) { - Expr *a = Expr::FromParam(param[0]); - Expr *b = Expr::FromParam(param[1]); - Expr *c = Expr::FromParam(param[2]); - Expr *d = Expr::FromParam(param[3]); - - Expr *two = Expr::FromConstant(2); - + ExprQuaternion q = (SS.GetEntity(normal))->NormalGetExprs(); // Convert the quaternion to our plane's normal vector. - n->x = two->Times(a->Times(c)); - n->x = (n->x)->Plus (two->Times(b->Times(d))); - n->y = two->Times(c->Times(d)); - n->y = (n->y)->Minus(two->Times(a->Times(b))); - n->z = a->Square(); - n->z = (n->z)->Minus(b->Square()); - n->z = (n->z)->Minus(c->Square()); - n->z = (n->z)->Plus (d->Square()); + *n = q.RotationN(); ExprVector p0 = SS.GetEntity(point[0])->PointGetExprs(); // The plane is n dot (p - p0) = 0, or @@ -489,7 +476,7 @@ void Entity::DrawOrGetDistance(int order) { int i, c = 20; Vector prev = u.ScaledBy(r).Plus(center); - for(i = 0; i <= c; i++) { + for(i = 1; i <= c; i++) { double phi = (2*PI*i)/c; Vector p = (u.ScaledBy(r*cos(phi))).Plus( v.ScaledBy(r*sin(phi))); diff --git a/expr.cpp b/expr.cpp index 2691c8a..01abd72 100644 --- a/expr.cpp +++ b/expr.cpp @@ -99,6 +99,24 @@ ExprVector ExprQuaternion::RotationV(void) { return v; } +ExprVector ExprQuaternion::RotationN(void) { + ExprVector n; + Expr *two = Expr::FromConstant(2); + + n.x = two->Times( w->Times(vy)); + n.x = (n.x)->Plus (two->Times(vx->Times(vz))); + + n.y = two->Times(vy->Times(vz)); + n.y = (n.y)->Minus(two->Times( w->Times(vx))); + + n.z = w->Square(); + n.z = (n.z)->Minus(vx->Square()); + n.z = (n.z)->Minus(vy->Square()); + n.z = (n.z)->Plus (vz->Square()); + + return n; +} + Expr *Expr::FromParam(hParam p) { Expr *r = AllocExpr(); r->op = PARAM; diff --git a/expr.h b/expr.h index 44adf2c..74c0614 100644 --- a/expr.h +++ b/expr.h @@ -139,6 +139,7 @@ public: ExprVector RotationU(void); ExprVector RotationV(void); + ExprVector RotationN(void); }; #endif diff --git a/glhelper.cpp b/glhelper.cpp index 8cbf29e..a3a3564 100644 --- a/glhelper.cpp +++ b/glhelper.cpp @@ -133,8 +133,9 @@ void glxFillPolygon(SPolygon *p) gluTessCallback(gt, GLU_TESS_COMBINE, (cf *)Combine); gluTessProperty(gt, GLU_TESS_WINDING_RULE, GLU_TESS_WINDING_ODD); - Vector normal = p->Normal(); + Vector normal = p->normal; glNormal3d(normal.x, normal.y, normal.z); + gluTessNormal(gt, normal.x, normal.y, normal.z); gluTessBeginPolygon(gt, NULL); for(i = 0; i < p->l.n; i++) { @@ -154,3 +155,34 @@ void glxFillPolygon(SPolygon *p) gluDeleteTess(gt); } +void glxMarkPolygonNormal(SPolygon *p) { + Vector tail = Vector::MakeFrom(0, 0, 0); + int i, j, cnt = 0; + // Choose some reasonable center point. + for(i = 0; i < p->l.n; i++) { + SContour *sc = &(p->l.elem[i]); + for(j = 0; j < (sc->l.n-1); j++) { + SPoint *sp = &(sc->l.elem[j]); + tail = tail.Plus(sp->p); + cnt++; + } + } + if(cnt == 0) return; + tail = tail.ScaledBy(1.0/cnt); + + Vector gn = SS.GW.projRight.Cross(SS.GW.projUp); + Vector tip = tail.Plus((p->normal).WithMagnitude(40/SS.GW.scale)); + Vector arrow = (p->normal).WithMagnitude(15/SS.GW.scale); + + glColor3d(1, 1, 0); + glBegin(GL_LINES); + glxVertex3v(tail); + glxVertex3v(tip); + glxVertex3v(tip); + glxVertex3v(tip.Minus(arrow.RotatedAbout(gn, 0.6))); + glxVertex3v(tip); + glxVertex3v(tip.Minus(arrow.RotatedAbout(gn, -0.6))); + glEnd(); + glEnable(GL_LIGHTING); +} + diff --git a/polygon.cpp b/polygon.cpp index f6fbf4e..54bba79 100644 --- a/polygon.cpp +++ b/polygon.cpp @@ -83,9 +83,33 @@ void SPolygon::MakeEdgesInto(SEdgeList *el) { } } -Vector SPolygon::Normal(void) { +Vector SPolygon::ComputeNormal(void) { if(l.n < 1) return Vector::MakeFrom(0, 0, 0); - return (l.elem[0]).Normal(); + return (l.elem[0]).ComputeNormal(); +} + +void SPolygon::FixContourDirections(void) { + // Outside curve looks counterclockwise, projected against our normal. + int i, j; + for(i = 0; i < l.n; i++) { + SContour *sc = &(l.elem[i]); + if(sc->l.n < 1) continue; + Vector pt = (sc->l.elem[0]).p; + + bool outer = true; + for(j = 0; j < l.n; j++) { + if(i == j) continue; + SContour *sct = &(l.elem[j]); + if(sct->ContainsPointProjdToNormal(normal, pt)) { + outer = !outer; + } + } + + bool clockwise = sc->IsClockwiseProjdToNormal(normal); + if(clockwise && outer || (!clockwise && !outer)) { + sc->Reverse(); + } + } } void SContour::MakeEdgesInto(SEdgeList *el) { @@ -98,23 +122,81 @@ void SContour::MakeEdgesInto(SEdgeList *el) { } } -Vector SContour::Normal(void) { - if(l.n < 3) return Vector::MakeFrom(0, 0, 0); +Vector SContour::ComputeNormal(void) { + Vector n = Vector::MakeFrom(0, 0, 0); - Vector u = (l.elem[0].p).Minus(l.elem[1].p); - - Vector v; - double dot = 2; - // Find the edge in the contour that's closest to perpendicular to the - // first edge, since that will give best numerical stability. - for(int i = 1; i < (l.n-1); i++) { - Vector vt = (l.elem[i].p).Minus(l.elem[i+1].p); - double dott = fabs(vt.Dot(u)/(u.Magnitude()*vt.Magnitude())); - if(dott < dot) { - dot = dott; - v = vt; + for(int i = 0; i < l.n - 2; i++) { + Vector u = (l.elem[i+1].p).Minus(l.elem[i+0].p).WithMagnitude(1); + Vector v = (l.elem[i+2].p).Minus(l.elem[i+1].p).WithMagnitude(1); + Vector nt = u.Cross(v); + if(nt.Magnitude() > n.Magnitude()) { + n = nt; } } - return (u.Cross(v)).WithMagnitude(1); + return n; } +bool SContour::IsClockwiseProjdToNormal(Vector n) { + if(n.Magnitude() < 0.01) oops(); + // An arbitrary 2d coordinate system that has n as its normal + Vector u = n.Normal(0); + Vector v = n.Normal(1); + + double area = 0; + for(int i = 0; i < (l.n - 1); i++) { + double u0 = (l.elem[i ].p).Dot(u); + double v0 = (l.elem[i ].p).Dot(v); + double u1 = (l.elem[i+1].p).Dot(u); + double v1 = (l.elem[i+1].p).Dot(v); + + area += ((v0 + v1)/2)*(u1 - u0); + } + + return (area < 0); +} + +bool SContour::ContainsPointProjdToNormal(Vector n, Vector p) { + Vector u = n.Normal(0); + Vector v = n.Normal(1); + + double up = p.Dot(u); + double vp = p.Dot(v); + + bool inside = false; + for(int i = 0; i < (l.n - 1); i++) { + double ua = (l.elem[i ].p).Dot(u); + double va = (l.elem[i ].p).Dot(v); + double ub = (l.elem[i+1].p).Dot(u); + double vb = (l.elem[i+1].p).Dot(v); + + // Write the parametric equation of the line, standardized so that + // t = 0 has smaller v than t = 1 + double u0, v0, du, dv; + + if(va < vb) { + u0 = ua; v0 = va; + du = (ub - ua); dv = (vb - va); + } else { + u0 = ub; v0 = vb; + du = (ua - ub); dv = (va - vb); + } + + if(dv == 0) continue; // intersects our horiz ray either 0 or 2 times + + double t = (vp - v0)/dv; + double ui = u0 + t*du; + if(ui > up && t >= 0 && t < 1) inside = !inside; + } + + return inside; +} + +void SContour::Reverse(void) { + int i; + for(i = 0; i < (l.n / 2); i++) { + int i2 = (l.n - 1) - i; + SPoint t = l.elem[i2]; + l.elem[i2] = l.elem[i]; + l.elem[i] = t; + } +} diff --git a/polygon.h b/polygon.h index 5043a2d..25d5651 100644 --- a/polygon.h +++ b/polygon.h @@ -57,17 +57,22 @@ public: SList l; void MakeEdgesInto(SEdgeList *el); - Vector Normal(void); + void Reverse(void); + Vector ComputeNormal(void); + bool IsClockwiseProjdToNormal(Vector n); + bool ContainsPointProjdToNormal(Vector n, Vector p); }; class SPolygon { public: SList l; + Vector normal; - Vector Normal(void); + Vector ComputeNormal(void); void AddEmptyContour(void); void AddPoint(Vector p); void MakeEdgesInto(SEdgeList *el); + void FixContourDirections(void); void Clear(void); }; diff --git a/sketch.cpp b/sketch.cpp index 15c8535..6a8708e 100644 --- a/sketch.cpp +++ b/sketch.cpp @@ -175,6 +175,7 @@ void Group::MakePolygons(void) { SEdge error; if(edges.AssemblePolygon(&poly, &error)) { polyError.yes = false; + poly.normal = poly.ComputeNormal(); faces.Add(&poly); } else { polyError.yes = true; @@ -199,8 +200,19 @@ void Group::MakePolygons(void) { // The bottom memset(&poly, 0, sizeof(poly)); if(!edges.AssemblePolygon(&poly, &error)) oops(); + Vector n = poly.ComputeNormal(); + if(translate.Dot(n) > 0) { + n = n.ScaledBy(-1); + } + poly.normal = n; + poly.FixContourDirections(); + poly.FixContourDirections(); faces.Add(&poly); + // Regenerate the edges, with the contour directions fixed up. + edges.l.Clear(); + poly.MakeEdgesInto(&edges); + // The sides int i; for(i = 0; i < edges.l.n; i++) { @@ -212,7 +224,9 @@ void Group::MakePolygons(void) { poly.AddPoint((edge->b).Plus(translate)); poly.AddPoint((edge->a).Plus(translate)); poly.AddPoint(edge->a); + poly.normal = ((edge->a).Minus(edge->b).Cross(n)).WithMagnitude(1); faces.Add(&poly); + edge->a = (edge->a).Plus(translate); edge->b = (edge->b).Plus(translate); } @@ -220,6 +234,7 @@ void Group::MakePolygons(void) { // The top memset(&poly, 0, sizeof(poly)); if(!edges.AssemblePolygon(&poly, &error)) oops(); + poly.normal = n.ScaledBy(-1); faces.Add(&poly); } } @@ -248,6 +263,15 @@ void Group::Draw(void) { glMaterialfv(GL_FRONT_AND_BACK, GL_AMBIENT_AND_DIFFUSE, vec); for(i = 0; i < faces.n; i++) { glxFillPolygon(&(faces.elem[i])); +#if 0 + // Debug stuff to show normals to the faces on-screen + glDisable(GL_LIGHTING); + glDisable(GL_DEPTH_TEST); + glxMarkPolygonNormal(&(faces.elem[i])); + glEnable(GL_LIGHTING); + glEnable(GL_DEPTH_TEST); +#endif + } glDisable(GL_LIGHTING); } diff --git a/sketch.h b/sketch.h index d8b42ce..015ce94 100644 --- a/sketch.h +++ b/sketch.h @@ -159,7 +159,13 @@ public: static const int POINT_XFRMD = 2010; static const int NORMAL_IN_3D = 3000; - static const int NORMAL_IN_2D = 3001; + static const int NORMAL_COPY = 3001; + // This is a normal that lies in a plane; so if the defining workplane + // has basis vectors uw, vw, nw, then + // n = (cos theta)*uw + (sin theta)*vw + // u = (sin theta)*uw - (cos theta)*vw + // v = nw + static const int NORMAL_IN_PLANE = 3002; static const int NORMAL_XFRMD = 3010; static const int WORKPLANE = 10000; diff --git a/solvespace.h b/solvespace.h index f7cc433..381bd11 100644 --- a/solvespace.h +++ b/solvespace.h @@ -67,6 +67,7 @@ void MemFree(void *p); // Utility functions that are provided in the platform-independent code. void glxVertex3v(Vector u); void glxFillPolygon(SPolygon *p); +void glxMarkPolygonNormal(SPolygon *p); void glxWriteText(char *str); void glxWriteTextRefCenter(char *str); void glxTranslatev(Vector u);