From d9317e79d641985c34fb994dcb646423017dcaf6 Mon Sep 17 00:00:00 2001 From: Akash Anand <akasha@iitk.ac.in> Date: Mon, 6 Aug 2007 18:02:36 +0000 Subject: [PATCH] *** empty log message *** --- contrib/FourierModel/FM_Reader.cpp | 2 +- .../ParaboloidProjectionSurface.cpp | 281 ++++++++++++++++++ .../ParaboloidProjectionSurface.h | 69 +++++ .../FourierModel/PlaneProjectionSurface.cpp | 181 +++++++++++ contrib/FourierModel/PlaneProjectionSurface.h | 73 +++++ 5 files changed, 605 insertions(+), 1 deletion(-) create mode 100755 contrib/FourierModel/ParaboloidProjectionSurface.cpp create mode 100755 contrib/FourierModel/ParaboloidProjectionSurface.h create mode 100755 contrib/FourierModel/PlaneProjectionSurface.cpp create mode 100755 contrib/FourierModel/PlaneProjectionSurface.h diff --git a/contrib/FourierModel/FM_Reader.cpp b/contrib/FourierModel/FM_Reader.cpp index d86e5cd366..84740443b3 100644 --- a/contrib/FourierModel/FM_Reader.cpp +++ b/contrib/FourierModel/FM_Reader.cpp @@ -19,11 +19,11 @@ FM_Reader::FM_Reader(const char* fn) //for (unsigned int i = 0; i < _nPatches; i++) { unsigned int i = 0; while (!InputFile.eof()) { - std::cout << "here\n"; char psName[32]; InputFile >> psName; std::cout << psName; int psTag; + std::cout << psTag; InputFile >> psTag; double origin[3]; InputFile >> origin[0] >> origin[1] >> origin[2]; diff --git a/contrib/FourierModel/ParaboloidProjectionSurface.cpp b/contrib/FourierModel/ParaboloidProjectionSurface.cpp new file mode 100755 index 0000000000..e4d5519412 --- /dev/null +++ b/contrib/FourierModel/ParaboloidProjectionSurface.cpp @@ -0,0 +1,281 @@ +#include "ParaboloidProjectionSurface.h" + +ParaboloidProjectionSurface::ParaboloidProjectionSurface +(int tag) : ProjectionSurface(1.) +{ + SetTag(tag); + SetName(std::string("Paraboloid")); + + twoPi_ = 2 * M_PI; + K_[0] = K_[1] = 0.5; + + numParameters_ = 2; + + O_[0] = O_[1] = O_[2] = 0.; + + E0_[0] = 0.; E0_[1] = 0.; E0_[2] = 1.; + E1_[0] = 1.; E1_[1] = 0.; E1_[2] = 0.; + E2_[0] = 0.; E2_[1] = 1.; E2_[2] = 0.; + + scale_[0] = scale_[1] = scale_[2] = 1.; +} + +ParaboloidProjectionSurface::ParaboloidProjectionSurface +(int tag, double O[3], double E0[3], double E1[3], double scale[3], + double K[2]) : ProjectionSurface(1.) +{ + SetTag(tag); + SetName(std::string("Paraboloid")); + + twoPi_ = 2 * M_PI; + K_[0] = K[0]; K_[1] = K[1]; + + O_[0] = O[0]; O_[1] = O[1]; O_[2] = O[2]; + + E0_[0] = E0[0]; E0_[1] = E0[1]; E0_[2] = E0[2]; + E1_[0] = E1[0]; E1_[1] = E1[1]; E1_[2] = E1[2]; + + E2_[0] = E0_[1] * E1_[2] - E0_[2] * E1_[1]; + E2_[1] = E0_[2] * E1_[0] - E0_[0] * E1_[2]; + E2_[2] = E0_[0] * E1_[1] - E0_[1] * E1_[0]; + + scale_[0] = scale[0]; scale_[1] = scale[1]; scale_[2] = scale[2]; +} + +void ParaboloidProjectionSurface::F +(double u, double v, double &x, double &y, double &z) +{ + x = O_[0] + K_[0] * v * (scale_[1] * E1_[0] * cos(twoPi_ * (u - 0.5)) + + scale_[2] * E2_[0] * sin(twoPi_ * (u - 0.5))) + + K_[1] * v * v * scale_[0] * E0_[0]; + y = O_[1] + K_[0] * v * (scale_[1] * E1_[1] * cos(twoPi_ * (u - 0.5)) + + scale_[2] * E2_[1] * sin(twoPi_ * (u - 0.5))) + + K_[1] * v * v * scale_[0] * E0_[1]; + z = O_[2] + K_[0] * v * (scale_[1] * E1_[2] * cos(twoPi_ * (u - 0.5)) + + scale_[2] * E2_[2] * sin(twoPi_ * (u - 0.5))) + + K_[1] * v * v * scale_[0] * E0_[2]; +} + +bool ParaboloidProjectionSurface::Inverse +(double x, double y, double z, double &u,double &v) +{ + double R[3]; + R[0] = x - O_[0]; + R[1] = y - O_[1]; + R[2] = z - O_[2]; + + double RdotT = 0., RdotNcT = 0.; + for (int i=0;i<3;i++) { + RdotT += R[i] * E1_[i]; + RdotNcT += R[i] * E2_[i]; + } + + RdotT /= scale_[1]; + RdotNcT /= scale_[2]; + + u = atan2(RdotNcT,RdotT); + u /= twoPi_; + u += 0.5; + + v = sqrt(RdotT * RdotT + RdotNcT * RdotNcT); + v /= K_[0]; +} + +void ParaboloidProjectionSurface::Dfdu +(double u, double v, double &x, double &y, double &z) +{ + x = twoPi_ * K_[0] * v * (- scale_[1] * E1_[0] * sin(twoPi_ * (u - 0.5)) + + scale_[2] * E2_[0] * cos(twoPi_ * (u - 0.5))); + y = twoPi_ * K_[0] * v * (- scale_[1] * E1_[1] * sin(twoPi_ * (u - 0.5)) + + scale_[2] * E2_[1] * cos(twoPi_ * (u - 0.5))); + z = twoPi_ * K_[0] * v * (- scale_[1] * E1_[2] * sin(twoPi_ * (u - 0.5)) + + scale_[2] * E2_[2] * cos(twoPi_ * (u - 0.5))); +} + +void ParaboloidProjectionSurface::Dfdv +(double u, double v, double &x, double &y, double &z) +{ + x = K_[0] * (scale_[1] * E1_[0] * cos(twoPi_ * (u - 0.5)) + + scale_[2] * E2_[0] * sin(twoPi_ * (u - 0.5))) + + 2. * K_[1] * v * scale_[0] * E0_[0]; + y = K_[0] * (scale_[1] * E1_[1] * cos(twoPi_ * (u - 0.5)) + + scale_[2] * E2_[1] * sin(twoPi_ * (u - 0.5))) + + 2. * K_[1] * v * scale_[0] * E0_[1]; + z = K_[0] * (scale_[1] * E1_[2] * cos(twoPi_ * (u - 0.5)) + + scale_[2] * E2_[2] * sin(twoPi_ * (u - 0.5))) + + 2. * K_[1] * v * scale_[0] * E0_[2]; +} + +void ParaboloidProjectionSurface::Dfdfdudu +(double u,double v, double &x, double &y, double &z) +{ + x = - twoPi_ * twoPi_ * K_[0] * v * + (scale_[1] * E1_[0] * cos(twoPi_ * (u - 0.5)) + + scale_[2] * E2_[0] * sin(twoPi_ * (u - 0.5))); + y = - twoPi_ * twoPi_ * K_[0] * v * + (scale_[1] * E1_[1] * cos(twoPi_ * (u - 0.5)) + + scale_[2] * E2_[1] * sin(twoPi_ * (u - 0.5))); + z = - twoPi_ * twoPi_ * K_[0] * v * + (scale_[1] * E1_[2] * cos(twoPi_ * (u - 0.5)) + + scale_[2] * E2_[2] * sin(twoPi_ * (u - 0.5))); +} + +void ParaboloidProjectionSurface::Dfdfdudv +(double u, double v, double &x, double &y, double &z) +{ + x = K_[0] * twoPi_ * (- scale_[1] * E1_[0] * sin(twoPi_ * (u - 0.5)) + + scale_[2] * E2_[0] * cos(twoPi_ * (u - 0.5))); + y = K_[0] * twoPi_ * (- scale_[1] * E1_[1] * sin(twoPi_ * (u - 0.5)) + + scale_[2] * E2_[1] * cos(twoPi_ * (u - 0.5))); + z = K_[0] * twoPi_ * (- scale_[1] * E1_[2] * sin(twoPi_ * (u - 0.5)) + + scale_[2] * E2_[2] * cos(twoPi_ * (u - 0.5))); +} + +void ParaboloidProjectionSurface::Dfdfdvdv +(double u, double v, double &x, double &y, double &z) +{ + x = 2. * K_[1] * scale_[0] * E0_[0]; + y = 2. * K_[1] * scale_[0] * E0_[1]; + z = 2. * K_[1] * scale_[0] * E0_[2]; +} + +void ParaboloidProjectionSurface::Dfdfdfdududu +(double u,double v,double &x,double &y,double &z) +{ + x = twoPi_ * twoPi_ * twoPi_ * K_[0] * v * + (scale_[1] * E1_[0] * sin(twoPi_ * (u - 0.5)) - + scale_[2] * E2_[0] * cos(twoPi_ * (u - 0.5))); + y = twoPi_ * twoPi_ * twoPi_ * K_[0] * v * + (scale_[1] * E1_[1] * sin(twoPi_ * (u - 0.5)) - + scale_[2] * E2_[1] * cos(twoPi_ * (u - 0.5))); + z = twoPi_ * twoPi_ * twoPi_ * K_[0] * v * + (scale_[1] * E1_[2] * sin(twoPi_ * (u - 0.5)) + + scale_[2] * E2_[2] * cos(twoPi_ * (u - 0.5))); +} + +void ParaboloidProjectionSurface::Dfdfdfdududv +(double u,double v,double &x,double &y,double &z) +{ + x = - twoPi_ * twoPi_ * K_[0] * + (scale_[1] * E1_[0] * cos(twoPi_ * (u - 0.5)) + + scale_[2] * E2_[0] * sin(twoPi_ * (u - 0.5))); + y = - twoPi_ * twoPi_ * K_[0] * + (scale_[1] * E1_[1] * cos(twoPi_ * (u - 0.5)) + + scale_[2] * E2_[1] * sin(twoPi_ * (u - 0.5))); + z = - twoPi_ * twoPi_ * K_[0] * + (scale_[1] * E1_[2] * cos(twoPi_ * (u - 0.5)) + + scale_[2] * E2_[2] * sin(twoPi_ * (u - 0.5))); +} + +void ParaboloidProjectionSurface::Dfdfdfdudvdv +(double u,double v,double &x,double &y,double &z) +{ + x = y = z = 0.; +} + +void ParaboloidProjectionSurface::Dfdfdfdvdvdv +(double u,double v,double &x,double &y,double &z) +{ + x = y = z = 0.; +} + +bool ParaboloidProjectionSurface::OrthoProjectionOnSurface +(double x, double y, double z, double& u, double& v) +{ + double R[3]; + R[0] = x - O_[0]; + R[1] = y - O_[1]; + R[2] = z - O_[2]; + + double RdotT = 0., RdotNcT = 0.; + for (int i=0;i<3;i++) { + RdotT += R[i] * E1_[i]; + RdotNcT += R[i] * E2_[i]; + } + + RdotT /= scale_[1]; + RdotNcT /= scale_[2]; + + u = atan2(RdotNcT,RdotT); + u /= twoPi_; + u += 0.5; + + double r1 = 0., r2 = 0.; + for (int i=0;i<3;i++) { + r1 += R[i] * (E1_[i] * scale_[1] * cos(twoPi_ * (u - 0.5)) + + E2_[i] * scale_[2] * sin(twoPi_ * (u - 0.5))); + r2 += R[i] * scale_[0] * E0_[i]; + } + + double A = + scale_[1] * scale_[1] * cos(twoPi_ * (u - 0.5)) * cos(twoPi_ * (u - 0.5)) + + scale_[2] * scale_[2] * sin(twoPi_ * (u - 0.5)) * sin(twoPi_ * (u - 0.5)); + double B = scale_[0] * scale_[0]; + + double a = 2 * K_[1] * K_[1] * B * B; + double b = K_[0] * K_[0] * A * A - 2 * K_[1] * r2 * B; + double c = - K_[0] * r1 * A; + + std::vector<double> root = SolveCubic(a,b,c); + + if (root.size()) { + /* + double xP,yP,zP; + double minDist = 1.e12; + for (int i=0;i<root.size();i++) { + F(u,root[i],xP,yP,zP); + double dist = sqrt((x-xP)*(x-xP)+(y-yP)*(y-yP)+(z-zP)*(z-zP)); + if (dist < minDist) { + minDist = dist; + v = root[i]; + } + } + */ + for (int i=0;i<root.size();i++) { + if (root[i] >= 0.) { + v = root[i]; + break; + } + } + double tol =1.e-4; + if ((u > - tol) && (u < 1. + tol) && (v > - tol) && (v < 1. + tol)) + return true; + else + return false; + } + else + return false; +} + +void ParaboloidProjectionSurface:: +SetParameter(int i, double x) +{ + switch (i) { + case 0: + K_[0] = x; + case 1: + K_[1] = x; + } +} + +double ParaboloidProjectionSurface:: +GetParameter(int i) +{ + switch (i) { + case 0: + return K_[0]; + case 1: + return K_[1]; + } +} + +std::string ParaboloidProjectionSurface:: +GetLabel(int i) +{ + switch (i) { + case 0: + return std::string("Parabola Scale 1"); + case 1: + return std::string("Parabola Scale 2"); + } +} diff --git a/contrib/FourierModel/ParaboloidProjectionSurface.h b/contrib/FourierModel/ParaboloidProjectionSurface.h new file mode 100755 index 0000000000..0d0a5f58ae --- /dev/null +++ b/contrib/FourierModel/ParaboloidProjectionSurface.h @@ -0,0 +1,69 @@ +#ifndef _PARABOLOID_PROJECTION_SURFACE_H_ +#define _PARABOLOID_PROJECTION_SURFACE_H_ + +#include <cmath> +#include <vector> +#include "Utils.h" +#include "ProjectionSurface.h" + +class ParaboloidProjectionSurface : public ProjectionSurface { + private: + double twoPi_; + double K_[2]; + public: + ParaboloidProjectionSurface + (int tag); + ParaboloidProjectionSurface + (int tag, double O[3], double E0[3], double E1[3], double scale[3], + double K[2]); + ParaboloidProjectionSurface(ParaboloidProjectionSurface *ps) + : ProjectionSurface(ps) + { + twoPi_ = ps->twoPi_; + K_[0] = ps->K_[0]; + K_[1] = ps->K_[1]; + } + + virtual ~ParaboloidProjectionSurface + () {} + + virtual ProjectionSurface *clone() + { + return new ParaboloidProjectionSurface(this); + } + + // Abstract methods of ProjectionSurface + + virtual void F + (double u, double v, double &x, double &y, double &z); + virtual bool Inverse + (double x,double y,double z,double &u,double &v); + virtual void Dfdu + (double u, double v, double &x, double &y, double &z); + virtual void Dfdv + (double u, double v, double &x, double &y, double &z); + virtual void Dfdfdudu + (double u,double v,double &x,double &y,double &z); + virtual void Dfdfdudv + (double u,double v,double &x,double &y,double &z); + virtual void Dfdfdvdv + (double u,double v,double &x,double &y,double &z); + virtual void Dfdfdfdududu + (double u,double v,double &x,double &y,double &z); + virtual void Dfdfdfdududv + (double u,double v,double &x,double &y,double &z); + virtual void Dfdfdfdudvdv + (double u,double v,double &x,double &y,double &z); + virtual void Dfdfdfdvdvdv + (double u,double v,double &x,double &y,double &z); + virtual bool OrthoProjectionOnSurface + (double x, double y, double z, double &u,double &v); + virtual void SetParameter + (int i, double x); + virtual double GetParameter + (int i); + virtual std::string GetLabel + (int i); +}; + +#endif diff --git a/contrib/FourierModel/PlaneProjectionSurface.cpp b/contrib/FourierModel/PlaneProjectionSurface.cpp new file mode 100755 index 0000000000..675967dab8 --- /dev/null +++ b/contrib/FourierModel/PlaneProjectionSurface.cpp @@ -0,0 +1,181 @@ +#include "PlaneProjectionSurface.h" + +PlaneProjectionSurface::PlaneProjectionSurface +(int tag) : ProjectionSurface(1.) +{ + SetTag(tag); + SetName(std::string("Plane")); + + numParameters_ = 0; + + O_[0] = O_[1] = O_[2] = 0.; + + E0_[0] = 1.; E0_[1] = 0.; E0_[2] = 0.; + E1_[0] = 0.; E1_[1] = 1.; E1_[2] = 0.; + E2_[0] = 0.; E2_[1] = 0.; E2_[2] = 1.; + + scale_[0] = scale_[1] = scale_[2] = 1.; +} + +PlaneProjectionSurface::PlaneProjectionSurface +(int tag, double O[3], double E0[3], double E1[3], double scale[3]) + : ProjectionSurface(1.) +{ + SetTag(tag); + SetName(std::string("Plane")); + + O_[0] = O[0]; O_[1] = O[1]; O_[2] = O[2]; + + E0_[0] = E0[0]; E0_[1] = E0[1]; E0_[2] = E0[2]; + E1_[0] = E1[0]; E1_[1] = E1[1]; E1_[2] = E1[2]; + + E2_[0] = E0_[1] * E1_[2] - E0_[2] * E1_[1]; + E2_[1] = E0_[2] * E1_[0] - E0_[0] * E1_[2]; + E2_[2] = E0_[0] * E1_[1] - E0_[1] * E1_[0]; + + scale_[0] = scale[0]; scale_[1] = scale[1]; scale_[2] = scale[2]; +} + +void PlaneProjectionSurface:: +F(double u, double v, double &x, double &y, double &z) +{ + x = O_[0] + u * scale_[1] * E1_[0] + v * scale_[2] * E2_[0]; + y = O_[1] + u * scale_[1] * E1_[1] + v * scale_[2] * E2_[1]; + z = O_[2] + u * scale_[1] * E1_[2] + v * scale_[2] * E2_[2]; +} + +bool PlaneProjectionSurface:: +Inverse(double x, double y, double z, double &u,double &v) +{ + u = (x - O_[0]) * E1_[0] + (y - O_[1]) * E1_[1] + (z - O_[2]) * E1_[2]; + v = (x - O_[0]) * E2_[0] + (y - O_[1]) * E2_[1] + (z - O_[2]) * E2_[2]; + + u /= scale_[1]; + v /= scale_[2]; + + double tol = 1.e-12; + + if ((u > - tol) && (u < 1. + tol) && (v > - tol) && (v < 1. + tol)) + return true; + else + return false; +} + +void PlaneProjectionSurface:: +Dfdu(double u, double v, double &x, double &y, double &z) +{ + x = scale_[1] * E1_[0]; + y = scale_[1] * E1_[1]; + z = scale_[1] * E1_[2]; +} + +void PlaneProjectionSurface:: +Dfdv(double u, double v, double &x, double &y, double &z) +{ + x = scale_[2] * E2_[0]; + y = scale_[2] * E2_[1]; + z = scale_[2] * E2_[2]; +} + +void PlaneProjectionSurface:: +Dfdfdudu(double u,double v, double &x, double &y, double &z) +{ + x = y = z = 0.; +} + +void PlaneProjectionSurface:: +Dfdfdudv(double u, double v, double &x, double &y, double &z) +{ + x = y = z = 0.; +} + +void PlaneProjectionSurface:: +Dfdfdvdv(double u, double v, double &x, double &y, double &z) +{ + x = y = z = 0.; +} + +void PlaneProjectionSurface:: +Dfdfdfdududu(double u,double v,double &x,double &y,double &z) +{ + x = y = z = 0.; +} + +void PlaneProjectionSurface:: +Dfdfdfdududv(double u,double v,double &x,double &y,double &z) +{ + x = y = z = 0.; +} + +void PlaneProjectionSurface:: +Dfdfdfdudvdv(double u,double v,double &x,double &y,double &z) +{ + x = y = z = 0.; +} + +void PlaneProjectionSurface:: +Dfdfdfdvdvdv(double u,double v,double &x,double &y,double &z) +{ + x = y = z = 0.; +} + +void PlaneProjectionSurface:: +GetNormal(double u, double v, double &x, double &y, double &z) +{ + x = E0_[0]; + y = E0_[1]; + z = E0_[2]; +} + +void PlaneProjectionSurface:: +Dndu(double u, double v, double &x, double &y, double &z) +{ + x = y = z = 0.; +} + +void PlaneProjectionSurface:: +Dndv(double u, double v, double &x, double &y, double &z) +{ + x = y = z = 0.; +} + +void PlaneProjectionSurface:: +Dndndudu(double u, double v, double &x, double &y, double &z) +{ + x = y = z = 0.; +} + +void PlaneProjectionSurface:: +Dndndudv(double u, double v, double &x, double &y, double &z) +{ + x = y = z = 0.; +} + +void PlaneProjectionSurface:: +Dndndvdv(double u, double v, double &x, double &y, double &z) +{ + x = y = z = 0.; +} + +bool PlaneProjectionSurface:: +OrthoProjectionOnSurface(double x, double y, double z, double &u,double &v) +{ + return Inverse(x,y,z,u,v); +} + +void PlaneProjectionSurface:: +SetParameter(int i, double x) +{ +} + +double PlaneProjectionSurface:: +GetParameter(int i) +{ + return 0.; +} + +std::string PlaneProjectionSurface:: +GetLabel(int i) +{ + return std::string(" "); +} diff --git a/contrib/FourierModel/PlaneProjectionSurface.h b/contrib/FourierModel/PlaneProjectionSurface.h new file mode 100755 index 0000000000..b592178996 --- /dev/null +++ b/contrib/FourierModel/PlaneProjectionSurface.h @@ -0,0 +1,73 @@ +#ifndef _PLANE_PROJECTION_SURFACE_H_ +#define _PLANE_PROJECTION_SURFACE_H_ + +#include <cmath> +#include "ProjectionSurface.h" + +class PlaneProjectionSurface : public ProjectionSurface { + public: + PlaneProjectionSurface + (int tag); + PlaneProjectionSurface + (int tag, double O[3], double E0[3], double E1[3], double scale[3]); + PlaneProjectionSurface(PlaneProjectionSurface *ps) + : ProjectionSurface(ps) {} + + virtual ~PlaneProjectionSurface + () {} + + virtual ProjectionSurface *clone() + { + return new PlaneProjectionSurface(this); + } + + // Abstract methods of ProjectionSurface + + virtual void F + (double u, double v, double &x, double &y, double &z); + virtual bool Inverse + (double x,double y,double z,double &u,double &v); + virtual void Dfdu + (double u, double v, double &x, double &y, double &z); + virtual void Dfdv + (double u, double v, double &x, double &y, double &z); + virtual void Dfdfdudu + (double u,double v,double &x,double &y,double &z); + virtual void Dfdfdudv + (double u,double v,double &x,double &y,double &z); + virtual void Dfdfdvdv + (double u,double v,double &x,double &y,double &z); + virtual void Dfdfdfdududu + (double u,double v,double &x,double &y,double &z); + virtual void Dfdfdfdududv + (double u,double v,double &x,double &y,double &z); + virtual void Dfdfdfdudvdv + (double u,double v,double &x,double &y,double &z); + virtual void Dfdfdfdvdvdv + (double u,double v,double &x,double &y,double &z); + virtual bool OrthoProjectionOnSurface + (double x, double y, double z, double &u,double &v); + virtual void SetParameter + (int i, double x); + virtual double GetParameter + (int i); + virtual std::string GetLabel + (int i); + + // Redefinitions for PlaneProjectionSurface + + virtual void GetNormal + (double u, double v, double &x, double &y, double &z); + virtual void Dndu + (double u, double v, double &x, double &y, double &z); + virtual void Dndv + (double u, double v, double &x, double &y, double &z); + virtual void Dndndudu + (double u, double v, double &x, double &y, double &z); + virtual void Dndndudv + (double u, double v, double &x, double &y, double &z); + virtual void Dndndvdv + (double u, double v, double &x, double &y, double &z); +}; + +#endif -- GitLab