ellipsoid_ump.cpp

曲面加工程序源码,用C＋＋编写的加工曲线曲面的程序

#if !defined(__ellipsoid_ump_cpp)
#define __ellipsoid_ump_cpp

//Difinition of ellipsoid_ump class    2000.07.13

#include <math.h>
#pragma hdrstop

#include "ellipsoid_ump.h"

////Construcor and destructor
ellipsoid_ump::ellipsoid_ump(){

	cz=13.098538; //[mm]
	cy=13.098538; //[mm]
	cx=666.6666666; //[mm]

   dx=388.4; //[mm]
   dy=0.; //[mm]
   dz=-10.6459550319; //[mm]

   rx_design=-0.6; //rotation around x-axis(form)  [deg]
   rx_setup =0.;//ritation araound x-asis(setup)  [deg]

   h_plane=-50.;

   z_top=100.;//=38.86-5;
};

ellipsoid_ump::~ellipsoid_ump(){

};
//equations for calcuration
inline bool ellipsoid_ump::InMirrorArea(double x,double y){

   //return  InMirrorBody(x,y,PlaneZ(x,y));
   return true;
};

inline bool ellipsoid_ump::InMirrorBody(double x,double y,double z){

    double t=(rx_setup+rx_design)*M_PI/180.;

    double xo=x*cos(t)+z*sin(t)+dx;
    double yo=y+dy;
    double zo=-x*sin(t)+z*cos(t)+dz;

    return ((xo/cx)*(xo/cx)+(yo/cy)*(yo/cy)+(zo/cz)*(zo/cz)<=1);
};

inline int ellipsoid_ump::AreaJudgement(double x,double y, double &z){

   //儚乕僋偺僄儕傾傪敾掕偟丄僄儕傾ID偲偦偺揰偱偺Z嵗昗傪曉偡
   //AreaID=   0 :Mirror
   //          1 :Plane
   //          2 :ztop

   int judge;

   if(InMirrorArea(x,y)) {

       z=MirrorZ(x,y);
       judge=0;

   }else{

       z=PlaneZ(x,y);
       judge=1;
   };


   if (z>=z_top) {

           z=z_top;
           judge=2;
   };

   return judge;
};

inline double ellipsoid_ump::PlaneZ(double x,double y){

    //double b=rx_setup*M_PI/180.;
    /*if (fabs(b)<1.e-6) {
       return h_plane;
    }else{
       return ((h_plane-y*cos(b))/sin(b));
    };      */

    return h_plane;
};

inline double ellipsoid_ump::MirrorZ(double x,double y){

   //2000.02.07 morishin
   double t=(rx_setup+rx_design)*M_PI/180.;

   double z=-(2*pow(cy,2)*(pow(cz,2)*dx*sin(t) + cos(t)*(pow(cx,2)*dz - (pow(cx,2) - pow(cz,2))*x*sin(t))) +
      sqrt(-4*pow(cy,2)*(pow(cx,2)*pow(cos(t),2) + pow(cz,2)*pow(sin(t),2))*
         (-(pow(cx,2)*pow(cy,2)*pow(cz,2)) + pow(cy,2)*pow(cz,2)*pow(dx,2) +
           pow(cx,2)*pow(cy,2)*pow(dz,2) + pow(cx,2)*pow(cz,2)*pow(y,2) +
           2*pow(cy,2)*pow(cz,2)*dx*x*cos(t) + pow(cy,2)*pow(cz,2)*pow(x,2)*pow(cos(t),2) -
           2*pow(cx,2)*pow(cy,2)*dz*x*sin(t) + pow(cx,2)*pow(cy,2)*pow(x,2)*pow(sin(t),2)) +
        4*pow(cy,4)*pow(pow(cz,2)*dx*sin(t) +
           cos(t)*(pow(cx,2)*dz - (pow(cx,2) - pow(cz,2))*x*sin(t)),2)))/
   (2.*pow(cy,2)*(pow(cx,2)*pow(cos(t),2) + pow(cz,2)*pow(sin(t),2)));

   return z;
};
inline double ellipsoid_ump::MirrorgradX(double x,double y){

   //2000.02.07 morishin
   double t=(rx_setup+rx_design)*M_PI/180.;

   double z=-(-(pow(cy,2)*(pow(cx,2) - pow(cz,2))*sin(2*t)) -
      (2*sqrt(2)*pow(cx,2)*pow(cy,4)*pow(cz,2)*(x + dx*cos(t) - dz*sin(t)))/
       sqrt(pow(cx,2)*pow(cy,2)*pow(cz,2)*
         (pow(cx,2)*pow(cy,2) + pow(cy,2)*pow(cz,2) - pow(cy,2)*pow(dx,2) -
           pow(cy,2)*pow(dz,2) - 2*pow(cy,2)*pow(x,2) - pow(cx,2)*pow(y,2) -
           pow(cz,2)*pow(y,2) - 4*pow(cy,2)*dx*x*cos(t) +
           (pow(cx,2)*pow(cy,2) - pow(cy,2)*pow(cz,2) - pow(cy,2)*pow(dx,2) +
              pow(cy,2)*pow(dz,2) - pow(cx,2)*pow(y,2) + pow(cz,2)*pow(y,2))*cos(2*t) +
           4*pow(cy,2)*dz*x*sin(t) + 2*pow(cy,2)*dx*dz*sin(2*t))))/
   (2.*pow(cy,2)*(pow(cx,2)*pow(cos(t),2) + pow(cz,2)*pow(sin(t),2)));

   return z;

};
inline double ellipsoid_ump::MirrorgradY(double x,double y){

   //2000.02.07 morishin
   double t=(rx_setup+rx_design)*M_PI/180.;

   double z=(2*pow(cx,2)*pow(cz,2)*y)/
   sqrt(-4*pow(cy,2)*(pow(cx,2)*pow(cos(t),2) + pow(cz,2)*pow(sin(t),2))*
      (-(pow(cx,2)*pow(cy,2)*pow(cz,2)) + pow(cy,2)*pow(cz,2)*pow(dx,2) +
        pow(cx,2)*pow(cy,2)*pow(dz,2) + pow(cx,2)*pow(cz,2)*pow(y,2) +
        2*pow(cy,2)*pow(cz,2)*dx*x*cos(t) + pow(cy,2)*pow(cz,2)*pow(x,2)*pow(cos(t),2) -
        2*pow(cx,2)*pow(cy,2)*dz*x*sin(t) + pow(cx,2)*pow(cy,2)*pow(x,2)*pow(sin(t),2)) +
     4*pow(cy,4)*pow(pow(cz,2)*dx*sin(t) + cos(t)*(pow(cx,2)*dz - (pow(cx,2) - pow(cz,2))*x*sin(t)),
       2));

   return z;

};
inline double ellipsoid_ump::MirrorgradXX(double x,double y){
   //2000.02.07 morishin
   double t=(rx_setup+rx_design)*M_PI/180.;

   double z=(2*sqrt(2)*pow(cx,4)*pow(cy,4)*pow(cz,4)*(pow(cy,2) - pow(y,2)))/
   pow(pow(cx,2)*pow(cy,2)*pow(cz,2)*
     (pow(cx,2)*pow(cy,2) + pow(cy,2)*pow(cz,2) - pow(cy,2)*pow(dx,2) - pow(cy,2)*pow(dz,2) -
       2*pow(cy,2)*pow(x,2) - pow(cx,2)*pow(y,2) - pow(cz,2)*pow(y,2) - 4*pow(cy,2)*dx*x*cos(t) +
       (pow(cx,2)*pow(cy,2) - pow(cy,2)*pow(cz,2) - pow(cy,2)*pow(dx,2) + pow(cy,2)*pow(dz,2) -
          pow(cx,2)*pow(y,2) + pow(cz,2)*pow(y,2))*cos(2*t) + 4*pow(cy,2)*dz*x*sin(t) +
       2*pow(cy,2)*dx*dz*sin(2*t)),1.5);

   return z;

};
inline double ellipsoid_ump::MirrorgradYY(double x,double y){
   //2000.02.07 morishin
   double t=(rx_setup+rx_design)*M_PI/180.;

   double z=(2*sqrt(2)*pow(cx,4)*pow(cy,4)*pow(cz,4)*y*(x + dx*cos(t) - dz*sin(t)))/
   pow(pow(cx,2)*pow(cy,2)*pow(cz,2)*
     (pow(cx,2)*pow(cy,2) + pow(cy,2)*pow(cz,2) - pow(cy,2)*pow(dx,2) - pow(cy,2)*pow(dz,2) -
       2*pow(cy,2)*pow(x,2) - pow(cx,2)*pow(y,2) - pow(cz,2)*pow(y,2) - 4*pow(cy,2)*dx*x*cos(t) +
       (pow(cx,2)*pow(cy,2) - pow(cy,2)*pow(cz,2) - pow(cy,2)*pow(dx,2) + pow(cy,2)*pow(dz,2) -
          pow(cx,2)*pow(y,2) + pow(cz,2)*pow(y,2))*cos(2*t) + 4*pow(cy,2)*dz*x*sin(t) +
       2*pow(cy,2)*dx*dz*sin(2*t)),1.5);
   return z;

};

inline double ellipsoid_ump::GetZ( double x, double y){

   double z;

   AreaJudgement(x,y,z);

   return z;

   //z=0 ,when x=y=0.

};

inline vector ellipsoid_ump::GetPosition( double x, double y){
	return ( vector( x, y, GetZ(x,y) ) );
};

inline double ellipsoid_ump::gradX(double x, double y){

   double z;
   double value;

   switch (AreaJudgement(x,y,z)){

       case 0:           //Mirror
           value = MirrorgradX(x,y);
           break;

       case 1:           //Plane
           value = 0;
           break ;

       case 2:           //ztop
           value = 0;
           break ;

       default:
           value = 0;

   };

   return value;

};

inline double ellipsoid_ump::gradY(double x, double y){

   double z;
   double value;

   double b;

   switch (AreaJudgement(x,y,z)){

       case 0:
           value = MirrorgradY(x,y);
           break;

       case 1:
           b=rx_setup*M_PI/180.;

           if (tan(b)!=0){
               value= -1./tan(b);     //PlanegradY
           }else{
               value= 0;
           };

           break ;

       case 2:
           value = 0;   //ztop
           break ;

       default:
           value = 0;

   };

   return value;

};

inline double ellipsoid_ump::gradXX(double x, double y){

   double z;
   double value;

   switch (AreaJudgement(x,y,z)){

       case 0:           //Mirror
           value = MirrorgradXX(x,y);
           break;

       case 1:           //Plane
           value = 0;
           break ;

       case 2:           //ztop
           value = 0;
           break ;

       default:
           value = 0;

   };

   return value;

};

inline double ellipsoid_ump::gradYY(double x, double y){


   double z;
   double value;

   switch (AreaJudgement(x,y,z)){

       case 0:           //Mirror
           value = MirrorgradYY(x,y);
           break;

       case 1:           //Plane
           value = 0;
           break ;

       case 2:           //ztop
           value = 0;
           break ;

       default:
           value = 0;

   };

   return value;
};

// This is the normal vector of surface on the point(x,y).
// This is toward the inside of ellipsoid_ump.
// This.length is 1.
vector ellipsoid_ump::NormalVector ( double x, double y){
	return ( vector( -gradX(x,y), -gradY(x,y), 1).normalize() );  //OK 980520
};
vector ellipsoid_ump::NormalVector ( vector position){
	return ( NormalVector( position.getX(), position.getY() ) );
};

////GradientVector////////////////////////////////////////////
// This is one of tangential vactors on the tangential plane.
// This is toward direction which has the largest gradient.
// This length is 1.
vector ellipsoid_ump::GradientVector(double x, double y){
	double gx = gradX( x,y);
	double gy = gradY( x,y);

	return ( vector( gx, gy, gx*gx+gy*gy ).normalize() );
};
vector ellipsoid_ump::GradientVector( vector position){
	return( GradientVector( position.getX(), position.getY() ) );
};

////NoGradientVector///////////////////////////////////////////
// This is one of tangential vactors on the tangential plane.
// This is toward direction which has no gradient.
//// This length will be 1,
//// and this Z-coordinate will be Zero.
vector ellipsoid_ump::NoGradientVector(double x,double y){
	return( NormalVector(x,y)%GradientVector(x,y) );
};
vector ellipsoid_ump::NoGradientVector(vector position){
	return( NormalVector(position)%GradientVector(position) );
};


int ellipsoid_ump::WriteFaceName(ofstream& fout){

	fout << "(Face name is : ellipsoid_ump);" << endl;
	return(1);
};
int ellipsoid_ump::WriteParameters(ofstream& fout){

	fout << "(cx[mm] is :" << cx <<");" << endl;
	fout << "(cy[mm] is :" << cy <<");" << endl;
	fout << "(cz[mm] is :" << cz <<");" << endl;

	fout << "(dx[mm] is :" << dx <<");" << endl;
	fout << "(dy[mm] is :" << dy <<");" << endl;
	fout << "(dz[mm] is :" << dz <<");" << endl;

	fout << "(rx_design is[deg] :" << rx_design <<");" << endl;
	fout << "(rx_setup[deg] is :" << rx_setup <<");" << endl;

	fout << "(h_plane[mm] is :" << h_plane <<");" << endl;
	fout << "(z_top[mm] is :" << z_top <<");" << endl;

	return(1);
};


void ellipsoid_ump::SetParameter(TStringList *str)
{
   for (int i=0; i<str->Count; i++){

         AnsiString temp;

         temp = "cx";
         if (str->Names[i]==temp){
              cx = str->Values[str->Names[i]].ToDouble();
         };
         temp = "cy";
         if (str->Names[i]==temp){
              cy = str->Values[str->Names[i]].ToDouble();
         };
         temp = "cz";
         if (str->Names[i]==temp){
              cz = str->Values[str->Names[i]].ToDouble();
         };

         temp = "dx";
         if (str->Names[i]==temp){
              dx = str->Values[str->Names[i]].ToDouble();
         };
         temp = "dy";
         if (str->Names[i]==temp){
              dy = str->Values[str->Names[i]].ToDouble();
         };
         temp = "dz";
         if (str->Names[i]==temp){
              dz = str->Values[str->Names[i]].ToDouble();
         };

         temp = "rx_design";
         if (str->Names[i]==temp){
              rx_design = str->Values[str->Names[i]].ToDouble();
         };
         temp = "rx_setup";
         if (str->Names[i]==temp){
              rx_setup = str->Values[str->Names[i]].ToDouble();
         };

         temp = "h_plane";
         if (str->Names[i]==temp){
              h_plane = str->Values[str->Names[i]].ToDouble();
         };
         temp = "z_top";
         if (str->Names[i]==temp){
              z_top = str->Values[str->Names[i]].ToDouble();
         };

         //isConcave
   };
   face::SetFaceParameter(str);
};

void ellipsoid_ump::GetParameter(TStringList *str)
{

   str->Add(AnsiString("---ellipsoid_ump parameter(s)---"));
   str->Add(AnsiString("FaceType=ellipsoid_ump"));

   str->Add(AnsiString("---(x/cx)^2+(y/cy)^2+(z/cz)^2=1---"));

   str->Add(AnsiString("cx=")+ cx);
   str->Add(AnsiString("cy=")+ cy);
   str->Add(AnsiString("cz=")+ cz);

   str->Add(AnsiString("dx=")+ dx);
   str->Add(AnsiString("dy=")+ dy);
   str->Add(AnsiString("dz=")+ dz);

   str->Add(AnsiString("rx_design=")+ rx_design);
   str->Add(AnsiString("rx_setup=")+ rx_setup);

   str->Add(AnsiString("h_plane=")+ h_plane);
   str->Add(AnsiString("z_top=")+ z_top);

   face::GetFaceParameter(str);
};

#endif