C++ SplineSurface类代码示例

int main(int argc, char* argv[] )
{

    ALWAYS_ERROR_IF(argc != 7, "Usage: " << argv[0]
		    << " curve1infile curve2infile surfaceoutfile point_x point_y point_z" << endl);

    // Open input curve 1 file
    ifstream is1(argv[1]);
    ALWAYS_ERROR_IF(is1.bad(), "Bad or no curve 1 input filename");

    // Open input curve 2 file
    ifstream is2(argv[2]);
    ALWAYS_ERROR_IF(is2.bad(), "Bad or no curve 2 input filename");

    // Open output surface file
    ofstream os(argv[3]);
    ALWAYS_ERROR_IF(os.bad(), "Bad output filename");

    Point pt(atof(argv[4]), atof(argv[5]), atof(argv[6]));

    // Read curves from file
    SplineCurve curv1, curv2;
    ObjectHeader head;
    is1 >> head >> curv1;
    is2 >> head >> curv2;

    SplineSurface* surf = SweepSurfaceCreator::linearSweptSurface(curv1, curv2, pt);

    surf->writeStandardHeader(os);
    surf->write(os);
}

//===========================================================================
SplineSurface* Torus::createNonRationalSpline(double eps) const
//===========================================================================
{
  // First fetch the first circular boundary curve in the minor
  // direction
  shared_ptr<Circle> circ = getMinorCircle(domain_.vmin());

  // Feth non-rational spline approximation
  shared_ptr<SplineCurve> crv(circ->createNonRationalSpline(0.5*eps));

  // Rotate this circle the valid angle around the main axis to 
  // create the non-rational spline surface
  // Note that the result will be rational if the tolerance is equal to zero
  int status;
  SISLCurve *qc = Curve2SISL(*crv);
  double *point = const_cast<double*>(location_.begin());
  double *axis =  const_cast<double*>(z_axis_.begin());
  SISLSurf *qs = NULL;
  s1302(qc, 0.5*eps, parbound_.umax()-parbound_.umin(), point, axis,
	&qs, &status);
  if (status < 0 || qs == NULL)
    return NULL;  // Approximation failed

  SplineSurface *surf = SISLSurf2Go(qs);
  surf->setParameterDomain(domain_.umin(), domain_.umax(),
			   domain_.vmin(), domain_.vmax());
  if (isSwapped())
    surf->swapParameterDirection();

  freeSurf(qs);
  return surf;
}

int main(int argc, char* argv[] )
{

    ALWAYS_ERROR_IF(argc != 10, "Usage: " << argv[0]
		    << " curveinfile surfaceoutfile angle point_x point_y point_z axis_x axis_y axis_z" << endl);

    // Open input curve file
    ifstream is(argv[1]);
    ALWAYS_ERROR_IF(is.bad(), "Bad or no curve input filename");

    // Open output surface file
    ofstream os(argv[2]);
    ALWAYS_ERROR_IF(os.bad(), "Bad output filename");

    double angle = atof(argv[3]);
    Point pt(atof(argv[4]), atof(argv[5]), atof(argv[6]));
    Point axis(atof(argv[7]), atof(argv[8]), atof(argv[9]));

    // Read curve from file
    SplineCurve curve;
    ObjectHeader head;
    is >> head >> curve;

    SplineSurface* surf = SweepSurfaceCreator::rotationalSweptSurface(curve, angle, pt, axis);

    surf->writeStandardHeader(os);
    surf->write(os);
}

//===========================================================================
SplineSurface*
LoftSurfaceCreator::loftSurfaceFromUnifiedCurves(vector<shared_ptr<SplineCurve> >::iterator
						 first_curve,
						 vector<double>::iterator first_param,
						 int nmb_crvs)
//===========================================================================
{

  SplineSurface* surf = loftNonrationalSurface(first_curve, first_param, nmb_crvs);
  if (first_curve[0]->rational())
    {
      int n = surf->numCoefs_u() * surf->numCoefs_v();
      int kdim = surf->dimension() + 1;
      bool all_positive = true;
      vector<double>::const_iterator it = surf->rcoefs_begin();
      it += (kdim - 1);
      for (int i = 0; i < n; ++i)
	if (it[kdim * i] <= 0.0)
	  {
	    all_positive = false;
	    break;
	  }
      if (!all_positive)
	{
	  delete surf;
	  surf = loftRationalSurface(first_curve, first_param, nmb_crvs);
	}
    }

  return surf;
}

//==========================================================================
void GeometryTools::splitSurfaceIntoPatches(const SplineSurface& sf,
			     vector<SplineSurface>& pat)
//==========================================================================
{
    SplineSurface orig = sf;
    orig.makeBernsteinKnotsU();
    orig.makeBernsteinKnotsV();

    int num_u = orig.numCoefs_u();
    int num_v = orig.numCoefs_v();
    int order_u = orig.order_u();
    int order_v = orig.order_v();
    int numpat_u = num_u / order_u;
    int numpat_v = num_v / order_v;

    pat.resize(numpat_u * numpat_v);
    typedef vector<double>::const_iterator const_iter;
    const_iter itu = orig.basis_u().begin();
    const_iter itv;
    for (int i = 0; i < numpat_u; ++i) {
	itv = orig.basis_v().begin();
	for (int j = 0; j < numpat_v; ++j) {
	    shared_ptr<SplineSurface>
		new_sf(orig.subSurface(*itu, *itv,
				       *(itu+order_u), *(itv+order_v)));
	    pat[numpat_u*j + i] = *new_sf;
	    itv += order_v;
	}
	itu += order_u;
    }

    return;
}

int main(int argc, char* argv[] )
{
    ALWAYS_ERROR_IF(argc < 3, "Usage: " << argv[0]
                    << " inputsurf inputpoints" << endl);


    // Open input surface file
    ifstream is(argv[1]);
    ALWAYS_ERROR_IF(is.bad(), "Bad or no input filename");

    // Read surface from file
    SplineSurface sf;
    is >> sf;

    // Get points
    ifstream pts(argv[2]);
    ALWAYS_ERROR_IF(pts.bad(), "Bad or no input filename");
    int n;
    pts >> n;
    vector<double> pt(n*2);
    for (int i = 0; i < n; ++i) {
        pts >> pt[2*i] >> pt[2*i+1];
    }

    std::vector<Point> p(3, Point(sf.dimension()));
    for (int i = 0; i < 10000; ++i) {
        for (int j = 0; j < n; ++j) {
            sf.point(p, pt[2*j], pt[2*j+1], 0);
        }
    }
//      cout << p[0] << p[1] << p[2] << (p[1] % p[2]);
}

//==========================================================================
void GeometryTools::findDominant(const SplineSurface& surface,
		  Vector3D& dominant_u, Vector3D& dominant_v)
//==========================================================================
{
    int nu = surface.numCoefs_u();
    int nv = surface.numCoefs_v();
    vector<double>::const_iterator start = surface.coefs_begin();
    Vector3D temp;
    // Dominant in u-direction
    dominant_u = Vector3D(0.0, 0.0, 0.0);
    for (int j = 0; j < nv; ++j) {
	for (int dd = 0; dd < 3; ++dd) {
	    temp[dd] = *(start + 3*(nu*j + (nu-1)) + dd)
		- *(start + 3*(nu*j) + dd);
	}
	dominant_u += temp;
    }
    // Dominant in v-direction
    dominant_v = Vector3D(0.0, 0.0, 0.0);
    for (int i = 0; i < nu; ++i) {
	for (int dd = 0; dd < 3; ++dd) {
	    temp[dd] = *(start + 3*(nu*(nv-1) + i) + dd)
		- *(start + 3*i + dd);
	}
	dominant_v += temp;
    }

    return;
}

//===========================================================================
BoundingBox Torus::boundingBox() const
//===========================================================================
{
    // A rather unefficient hack...
    SplineSurface* tmp = geometrySurface();
    BoundingBox box = tmp->boundingBox();
    delete tmp;
    return box;
}

int main()
{
    ObjectHeader header;
    SplineSurface surf;
    cin >> header >> surf;

    PointCloud<3> cloud(surf.coefs_begin(),
                        surf.numCoefs_u()*surf.numCoefs_v());
    cloud.writeStandardHeader(cout);
    cout << cloud;
}

//===========================================================================
shared_ptr<SplineSurface> SplineUtils::refineToBezier(const SplineSurface& spline_sf)
//===========================================================================
{
    shared_ptr<SplineSurface> bez_sf;

    const BsplineBasis& bas_u = spline_sf.basis_u();
    const BsplineBasis& bas_v = spline_sf.basis_v();
    const int order_u = bas_u.order();
    const int order_v = bas_v.order();

    // We extract the unique knots.
    vector<double> new_knots_u, new_knots_v;
    // vector<double> ref_knots_u, ref_knots_v;
    vector<double>::const_iterator iter = bas_u.begin();
    while (iter != bas_u.end() - order_u)
    {
	if (iter[0] != iter[1])
	{
	    int knot_mult = bas_u.knotMultiplicity(iter[0]);
	    int num_insert = order_u - knot_mult;
	    if (num_insert > 0)
	    {
		new_knots_u.insert(new_knots_u.end(), num_insert, iter[0]);
	    }
	}
	++iter;
    }

    iter = bas_v.begin();
    while (iter != bas_v.end() - order_v)
    {
	if (iter[0] != iter[1])
	{
	    int knot_mult = bas_v.knotMultiplicity(iter[0]);
	    int num_insert = order_v - knot_mult;
	    if (num_insert > 0)
	    {
		new_knots_v.insert(new_knots_v.end(), num_insert, iter[0]);
	    }
	}
	++iter;
    }

    bez_sf = insertKnots(spline_sf,
			 new_knots_u, new_knots_v);

    return bez_sf;
}

// ============================================================================
vector<CurveVec> SSurfTraceIsocontours(const SplineSurface& ss,
				       const vector<double>& isovals,
				       const double tol, 
				       bool include_3D_curves,
				       bool use_sisl_marching)
// ============================================================================
{
  assert(ss.dimension() == 1); // only intended to work for spline functions
  
  // Compute topology for each requested level-set.  We use SISL for this
  SISLSurf* sislsurf1D = GoSurf2SISL(ss, false);
  SISLSurf* sislsurf3D = use_sisl_marching ? make_sisl_3D(ss) : nullptr;
  
  // Defining function tracing out the level set for a specified isovalue
  const function<CurveVec(double)> comp_lset = [&] (double ival)
    {return compute_levelset(ss, sislsurf1D, sislsurf3D, ival, tol,
			     include_3D_curves, use_sisl_marching);};

  // Computing all level-set curves for all isovalues ("transforming" each
  // isovalue into its corresponding level-set)
  const vector<CurveVec> result = apply_transform(isovals, comp_lset);

  // Cleaning up after use of SISL objects
  freeSurf(sislsurf1D);
  if (sislsurf3D)
    freeSurf(sislsurf3D);

  // Returning result
  return result;
}

//==========================================================================
void create_bary_coord_system3D(const SplineSurface& surface,
				BaryCoordSystem3D& bc)
//==========================================================================
{
    BoundingBox box = surface.boundingBox();
    create_bary_coord_system3D(box, bc);
    return;
}

int main(int argc, char** argv)
{
    if (argc < 3) {
	cerr << "Usage: " << argv[0] << " u_res v_res" << endl;
	return 1;
    }
    int ures = atoi(argv[1]);
    int vres = atoi(argv[2]);

    ObjectHeader head;
    cin >> head;
    ASSERT(head.classType() == SplineSurface::classType());
    SplineSurface sf;
    cin >> sf;
    vector<double> points;
    vector<double> param_u;
    vector<double> param_v;
    sf.gridEvaluator(ures, vres, points, param_u, param_v);
    RectGrid grid(ures, vres, sf.dimension(), &points[0]);
    grid.writeStandardHeader(cout);
    grid.write(cout);
}

int main(int argc, char** argv)
{
    // Read the curve from file
    std::ifstream input(argv[1]);
    if (input.bad()) {
        std::cerr << "File error (no file or corrupt file specified)."
                  << std::endl;
        return 1;
    }
    ObjectHeader header;
    SplineSurface surface;
    input >> header >> surface;

    // Loop through parameter space
    const int samples = 50;
    double increment_u = (surface.endparam_u()
                          - surface.startparam_u()) / (samples-1);
    double increment_v = (surface.endparam_v()
                          - surface.startparam_v()) / (samples-1);
    Point result;
    double param_u = surface.startparam_u();
    int prec = (int)std::cout.precision(15);
    for (int i = 0; i < samples; ++i) {
        double param_v = surface.startparam_v();
        for (int j = 0; j < samples; ++j) {
            surface.point(result, param_u, param_v);
            std::cout << result[0] << "\t" << result[1] << "\t"
                      << result[2] << std::endl;
            param_v += increment_v;
        }
        std::cout << std::endl;
        param_u += increment_u;
    }
    std::cout.precision(prec);

    return 0;
}

//===========================================================================
shared_ptr<SplineSurface> SurfaceCreators::mergeRationalParts(const SplineSurface& nom_sf,
							      const SplineSurface& den_sf,
							      bool weights_in_first)
//===========================================================================
{
    ASSERT((!nom_sf.rational()) && (!den_sf.rational()));
    ASSERT(den_sf.dimension() == 1);

    int dim = nom_sf.dimension();

    // We first make sure they share spline space.
    vector<shared_ptr<SplineSurface> > sfs;
    sfs.push_back(shared_ptr<SplineSurface>(nom_sf.clone()));
    sfs.push_back(shared_ptr<SplineSurface>(den_sf.clone()));
    double knot_diff_tol = 1e-06;
    GeometryTools::unifySurfaceSplineSpace(sfs, knot_diff_tol);

    vector<double> rcoefs;
    vector<double>::const_iterator iter = sfs[0]->coefs_begin();
    vector<double>::const_iterator riter = sfs[1]->coefs_begin();
    int num_coefs = sfs[0]->numCoefs_u()*sfs[0]->numCoefs_v();
    for (int ki = 0; ki < num_coefs; ++ki) {
	for (int kj = 0; kj < dim; ++kj) {
	    if (weights_in_first) {
		rcoefs.push_back(iter[ki*dim+kj]);
	    } else {
		rcoefs.push_back(iter[ki*dim+kj]*riter[ki]);
	    }
	}
	rcoefs.push_back(riter[ki]);
    }

    shared_ptr<SplineSurface> rat_sf(new SplineSurface
				     (sfs[0]->numCoefs_u(), sfs[0]->numCoefs_v(),
				      sfs[0]->order_u(), sfs[0]->order_v(),
				      sfs[0]->basis_u().begin(), sfs[0]->basis_v().begin(),
				      rcoefs.begin(), dim, true));

    return rat_sf;
}