//=========================================================================== // GoTools Core - SINTEF Geometry Tools Core library, version 2.0.1 // // Copyright (C) 2000-2007, 2010 SINTEF ICT, Applied Mathematics, Norway. // // This program is free software; you can redistribute it and/or // modify it under the terms of the GNU General Public License // as published by the Free Software Foundation version 2 of the License. // // This program is distributed in the hope that it will be useful, // but WITHOUT ANY WARRANTY; without even the implied warranty of // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the // GNU General Public License for more details. // // You should have received a copy of the GNU General Public License // along with this program; if not, write to the Free Software // Foundation, Inc., // 59 Temple Place - Suite 330, // Boston, MA 02111-1307, USA. // // Contact information: E-mail: tor.dokken@sintef.no // SINTEF ICT, Department of Applied Mathematics, // P.O. Box 124 Blindern, // 0314 Oslo, Norway. // // Other licenses are also available for this software, notably licenses // for: // - Building commercial software. // - Building software whose source code you wish to keep private. //=========================================================================== 00015 #ifndef _SPLINESURFACE_H 00016 #define _SPLINESURFACE_H 00017 00018 00019 #include "GoTools/geometry/ParamSurface.h" 00020 #include "GoTools/geometry/BsplineBasis.h" 00021 #include "GoTools/geometry/RectDomain.h" 00022 #include "GoTools/utils/ScratchVect.h" 00023 #include "GoTools/utils/config.h" 00024 00025 namespace Go 00026 { 00027 00028 class Interpolator; 00029 class SplineCurve; 00030 class DirectionCone; 00031 00032 struct BasisPtsSf 00033 { 00034 double param[2]; // Parameter tripple in which the basis functions are evaluated 00035 int left_idx[2]; // Index of the first non-zero knot interval for all parameter directions 00036 // The value of all basis functions, size equal to (degree_u+1)*(degree_v+1)*(degree_w+1) 00037 std::vector< double > basisValues; 00038 00039 void preparePts(double u, double v, int idx_u, int idx_v, int size) 00040 { 00041 param[0] = u; 00042 param[1] = v; 00043 left_idx[0] = idx_u; 00044 left_idx[1] = idx_v; 00045 basisValues.resize(size); 00046 } 00047 }; 00048 00049 struct BasisDerivsSf 00050 { 00051 double param[2]; // Parameter double in which the basis functions are evaluated 00052 int left_idx[2]; // Index of the first non-zero knot interval for all parameter directions 00053 // The value of all basis functions, size equal to (degree_u+1)*(degree_v+1) 00054 std::vector< double > basisValues; 00055 // the derivative of all basis functions in u direction, same size as previous 00056 std::vector< double > basisDerivs_u; 00057 // the derivative of all basis functions in v direction, same size as previous 00058 std::vector< double > basisDerivs_v; 00059 00060 void prepareDerivs(double u, double v, int idx_u, int idx_v, int size) 00061 { 00062 param[0] = u; 00063 param[1] = v; 00064 left_idx[0] = idx_u; 00065 left_idx[1] = idx_v; 00066 basisValues.resize(size); 00067 basisDerivs_u.resize(size); 00068 basisDerivs_v.resize(size); 00069 } 00070 }; 00071 00072 struct BasisDerivsSf2 00073 { 00074 double param[2]; // Parameter double in which the basis functions are evaluated 00075 int left_idx[2]; // Index of the first non-zero knot interval for all parameter directions 00076 // The value of all basis functions, size equal to (degree_u+1)*(degree_v+1) 00077 std::vector< double > basisValues; 00078 00079 // the derivative of all basis functions in u direction, same size as previous 00080 std::vector< double > basisDerivs_u; 00081 // the derivative of all basis functions in v direction, same size as previous 00082 std::vector< double > basisDerivs_v; 00083 00084 // the second derivative of all basis functions twice in u direction, same size as previous 00085 std::vector< double > basisDerivs_uu; 00086 // the second derivative of all basis functions in u and v direction, same size as previous 00087 std::vector< double > basisDerivs_uv; 00088 // the second derivative of all basis functions twice in v direction, same size as previous 00089 std::vector< double > basisDerivs_vv; 00090 00091 void prepareDerivs(double u, double v, int idx_u, int idx_v, 00092 int size) 00093 { 00094 param[0] = u; 00095 param[1] = v; 00096 left_idx[0] = idx_u; 00097 left_idx[1] = idx_v; 00098 basisValues.resize(size); 00099 basisDerivs_u.resize(size); 00100 basisDerivs_v.resize(size); 00101 basisDerivs_uu.resize(size); 00102 basisDerivs_uv.resize(size); 00103 basisDerivs_vv.resize(size); 00104 } 00105 }; 00106 00130 class GO_API SplineSurface : public ParamSurface 00131 { 00132 public: 00135 SplineSurface() 00136 : dim_(-1), rational_(false) 00137 { 00138 } 00139 00163 template <typename RandomIterator1, 00164 typename RandomIterator2, 00165 typename RandomIterator3> 00166 SplineSurface(int number1, 00167 int number2, 00168 int order1, 00169 int order2, 00170 RandomIterator1 knot1start, 00171 RandomIterator2 knot2start, 00172 RandomIterator3 coefsstart, 00173 int dim, 00174 bool rational = false) 00175 : dim_(dim), rational_(rational), 00176 basis_u_(number1, order1, knot1start), 00177 basis_v_(number2, order2, knot2start) 00178 { 00179 if (rational) { 00180 int n = (dim+1)*number1*number2; 00181 rcoefs_.resize(n); 00182 std::copy(coefsstart, coefsstart + n, rcoefs_.begin()); 00183 coefs_.resize(dim*number1*number2); 00184 updateCoefsFromRcoefs(); 00185 } else { 00186 int n = dim*number1*number2; 00187 coefs_.resize(n); 00188 std::copy(coefsstart, coefsstart + n, coefs_.begin()); 00189 } 00190 } 00191 00208 template <typename RandomIterator> 00209 SplineSurface(const BsplineBasis& basis_u, 00210 const BsplineBasis& basis_v, 00211 RandomIterator coefsstart, 00212 int dim, 00213 bool rational = false) 00214 : dim_(dim), rational_(rational), 00215 basis_u_(basis_u), 00216 basis_v_(basis_v) 00217 { 00218 int number1 = basis_u.numCoefs(); 00219 int number2 = basis_v.numCoefs(); 00220 if (rational) { 00221 int n = (dim+1)*number1*number2; 00222 rcoefs_.resize(n); 00223 std::copy(coefsstart, coefsstart + n, rcoefs_.begin()); 00224 coefs_.resize(dim*number1*number2); 00225 updateCoefsFromRcoefs(); 00226 } else { 00227 int n = dim*number1*number2; 00228 coefs_.resize(n); 00229 std::copy(coefsstart, coefsstart + n, coefs_.begin()); 00230 } 00231 } 00232 00234 virtual ~SplineSurface(); 00235 00236 // inherited from Streamable 00237 virtual void read (std::istream& is); 00238 00239 // inherited from Streamable 00240 virtual void write (std::ostream& os) const; 00241 00242 // inherited from GeomObject 00243 virtual BoundingBox boundingBox() const; 00244 00245 // inherited from GeomObject 00246 virtual int dimension() const; 00247 00249 void swap(SplineSurface& other); 00250 00251 // inherited from GeomObject 00252 virtual ClassType instanceType() const; 00253 00254 // inherited from GeomObject 00255 static ClassType classType() 00256 { return Class_SplineSurface; } 00257 // #if ((_MSC_VER > 0) && (_MSC_VER < 1300)) 00258 // virtual GeomObject* clone() const 00259 // { return new SplineSurface(*this); } 00260 // #else 00261 virtual SplineSurface* clone() const 00262 { return new SplineSurface(*this); } 00263 // #endif 00264 00266 virtual SplineSurface* asSplineSurface() 00267 { 00268 return this; 00269 } 00270 00271 // inherited from ParamSurface 00272 // #if ((_MSC_VER > 0) && (_MSC_VER < 1300)) 00273 // virtual const Domain& parameterDomain() const; 00274 // #else 00275 virtual const RectDomain& parameterDomain() const; 00276 // #endif 00277 00278 // inherited from ParamSurface 00279 virtual RectDomain containingDomain() const; 00280 00282 virtual bool inDomain(double u, double v) const; 00283 00284 virtual Point closestInDomain(double u, double v) const; 00285 00286 // inherited from ParamSurface 00287 virtual CurveLoop outerBoundaryLoop(double degenerate_epsilon 00288 = DEFAULT_SPACE_EPSILON) const; 00289 // inherited from ParamSurface 00290 virtual std::vector<CurveLoop> allBoundaryLoops(double degenerate_epsilon 00291 = DEFAULT_SPACE_EPSILON) const; 00292 00293 // inherited from ParamSurface 00294 virtual void point(Point& pt, double upar, double vpar) const; 00295 00296 // Output: Partial derivatives up to order derivs (pts[0]=S(u,v), 00297 // pts[1]=dS/du=S_u, pts[2]=S_v, pts[3]=S_uu, pts[4]=S_uv, pts[5]=S_vv, ...) 00298 // inherited from ParamSurface 00299 virtual void point(std::vector<Point>& pts, 00300 double upar, double vpar, 00301 int derivs, 00302 bool u_from_right = true, 00303 bool v_from_right = true, 00304 double resolution = 1.0e-12) const; 00305 00308 virtual double startparam_u() const; 00309 00312 virtual double startparam_v() const; 00313 00316 virtual double endparam_u() const; 00317 00320 virtual double endparam_v() const; 00321 00322 // inherited from ParamSurface 00323 virtual void normal(Point& n, double upar, double vpar) const; 00324 00326 enum NormalConeMethod { 00327 SederbergMeyers = 0, 00328 SMCornersFirst = 1 00329 }; 00333 DirectionCone normalCone(NormalConeMethod method) const; 00334 00340 virtual DirectionCone normalCone() const; 00341 00342 // inherited from ParamSurface 00343 virtual DirectionCone tangentCone(bool pardir_is_u) const; 00344 00345 // Creates a composite box enclosing the surface. The composite box 00346 // consists of an inner and an edge box. The inner box is 00347 // made from the interior control points of the surface, while the 00348 // edge box is made from the boundary control points. 00349 // Inherited from ParamSurface 00350 virtual CompositeBox compositeBox() const; 00351 00353 SplineSurface* normal() const; 00354 00362 SplineSurface* normalSurface() const; 00363 00372 SplineSurface* derivSurface(int ider1, int ider2) const; 00373 00383 SplineSurface* subSurface(double from_upar, 00384 double from_vpar, 00385 double to_upar, 00386 double to_vpar, 00387 double fuzzy = 00388 DEFAULT_PARAMETER_EPSILON) const; 00389 00390 // inherited from ParamSurface 00391 virtual std::vector<boost::shared_ptr<ParamSurface> > 00392 subSurfaces(double from_upar, double from_vpar, 00393 double to_upar, double to_vpar, 00394 double fuzzy = DEFAULT_PARAMETER_EPSILON) const; 00395 00397 virtual SplineSurface* mirrorSurface(const Point& pos, const Point& norm) const; 00398 00399 // inherited from ParamSurface 00400 virtual void closestPoint(const Point& pt, 00401 double& clo_u, 00402 double& clo_v, 00403 Point& clo_pt, 00404 double& clo_dist, 00405 double epsilon, 00406 const RectDomain* domain_of_interest = NULL, 00407 double *seed = 0) const; 00408 00409 // inherited from ParamSurface 00410 virtual void closestBoundaryPoint(const Point& pt, 00411 double& clo_u, 00412 double& clo_v, 00413 Point& clo_pt, 00414 double& clo_dist, 00415 double epsilon, 00416 const RectDomain* rd = NULL, 00417 double *seed = 0) const; 00418 00442 void appendSurface(ParamSurface* sf, int join_dir, 00443 int continuity, double& dist, bool repar=true); 00444 00445 void appendSurface(SplineSurface* sf, int join_dir, int continuity, double& dist); 00446 00458 void appendSurface(ParamSurface* sf, int join_dir, bool repar=true); 00459 00460 // inherited from ParamSurface 00461 virtual void getBoundaryInfo(Point& pt1, Point& pt2, 00462 double epsilon, SplineCurve*& cv, 00463 SplineCurve*& crosscv, double knot_tol = 1e-05) const; 00464 00482 void getBoundaryInfo(double par1, double par2, 00483 int bdindex, SplineCurve*& cv, 00484 SplineCurve*& crosscv, double knot_tol = 1e-05) const; 00485 00512 void getBoundaryIdx(Point& pt1, Point& pt2, 00513 double epsilon, int& bdindex, 00514 double& par1, double& par2, double knot_tol = 1e-05) const; 00515 00517 // corresponding boundary. The Boundaries are indexed as in \getBoundaryIdx(). 00527 void getBoundaryIdx(Point& pt1, double epsilon, int& bdindex, 00528 double knot_tol = 1e-05) const; 00529 00530 // inherited from ParamSurface 00531 virtual void turnOrientation(); 00532 00533 // inherited from ParamSurface 00534 virtual void swapParameterDirection(); 00535 00536 // inherited from ParamSurface 00537 virtual void reverseParameterDirection(bool direction_is_u); 00538 00544 00545 00546 00547 00548 00549 00550 int boundaryIndex(Point& param_pt1, Point& param_pt2) const; 00551 00557 virtual double area(double tol) const; 00558 00561 const BsplineBasis& basis_u() const 00562 { return basis_u_; } 00563 00566 const BsplineBasis& basis_v() const 00567 { return basis_v_; } 00568 00571 BsplineBasis& basis_u() 00572 { return basis_u_; } 00573 00576 BsplineBasis& basis_v() 00577 { return basis_v_; } 00578 00583 const BsplineBasis& basis(int i) const 00584 { return (i==0) ? basis_u_ : basis_v_; } 00585 00588 int numCoefs_u() const 00589 { return basis_u_.numCoefs(); } 00590 00593 int numCoefs_v() const 00594 { return basis_v_.numCoefs(); } 00595 00598 int order_u() const 00599 { return basis_u_.order(); } 00600 00603 int order_v() const 00604 { return basis_v_.order(); } 00605 00608 bool rational() const 00609 { return rational_; } 00610 00615 std::vector<double>::iterator coefs_begin() 00616 { return coefs_.begin(); } 00617 00622 std::vector<double>::iterator coefs_end() 00623 { return coefs_.end(); } 00624 00629 std::vector<double>::const_iterator coefs_begin() const 00630 { return coefs_.begin(); } 00631 00636 std::vector<double>::const_iterator coefs_end() const 00637 { return coefs_.end(); } 00638 00643 std::vector<double>::iterator rcoefs_begin() 00644 { return rcoefs_.begin(); } 00645 00650 std::vector<double>::iterator rcoefs_end() 00651 { return rcoefs_.end(); } 00652 00657 std::vector<double>::const_iterator rcoefs_begin() const 00658 { return rcoefs_.begin(); } 00659 00664 std::vector<double>::const_iterator rcoefs_end() const 00665 { return rcoefs_.end(); } 00666 00671 std::vector<double>::iterator ctrl_begin() 00672 { return rational_ ? rcoefs_.begin() : coefs_.begin(); } 00673 00678 std::vector<double>::iterator ctrl_end() 00679 { return rational_ ? rcoefs_.end() : coefs_.end(); } 00680 00685 std::vector<double>::const_iterator ctrl_begin() const 00686 { return rational_ ? rcoefs_.begin() : coefs_.begin(); } 00687 00692 std::vector<double>::const_iterator ctrl_end() const 00693 { return rational_ ? rcoefs_.end() : coefs_.end(); } 00694 00697 void getWeights(std::vector<double>& weights) const; 00698 00699 // inherited from ParamSurface 00700 virtual bool isDegenerate(bool& b, bool& r, 00701 bool& t, bool& l, double tolerance) const; 00702 00704 virtual void getDegenerateCorners(std::vector<Point>& deg_corners, double tol) const; 00705 00711 void setParameterDomain(double u1, double u2, double v1, double v2); 00712 00715 void insertKnot_u(double apar); 00716 00720 void insertKnot_u(const std::vector<double>& new_knots); 00721 00724 void insertKnot_v(double apar); 00725 00729 void insertKnot_v(const std::vector<double>& new_knots); 00730 00733 void removeKnot_u(double tpar); 00734 00737 void removeKnot_v(double tpar); 00738 00741 void makeBernsteinKnotsU(); 00742 00745 void makeBernsteinKnotsV(); 00746 00748 void makeSurfaceKRegular(); 00749 00753 int numberOfPatches_u() const; 00754 00758 int numberOfPatches_v() const; 00759 00763 double knotSpan(int pardir, int iknot) const; 00764 00770 void raiseOrder(int raise_u, int raise_v); 00771 00780 SplineCurve* constParamCurve(double parameter, 00781 bool pardir_is_u) const; 00782 00796 void constParamCurve(double parameter, 00797 bool pardir_is_u, 00798 SplineCurve*& cv, 00799 SplineCurve*& crosscv) const; 00800 00813 void getConstParamCurves(const std::vector<double>& params_u, 00814 const std::vector<double>& params_v, 00815 std::vector<boost::shared_ptr<SplineCurve> >& curves_u, 00816 std::vector<boost::shared_ptr<SplineCurve> >& curves_v); 00817 00818 // inherited from ParamSurface 00819 virtual std::vector< boost::shared_ptr<ParamCurve> > 00820 constParamCurves(double parameter, bool pardir_is_u) const; 00821 00830 00831 00843 // @@sbr What about the orientation? Seems it follows that of the sf. 00844 SplineCurve* edgeCurve(int ccw_edge_number) const; 00845 00860 void interpolate(Interpolator& interpolator1, 00861 Interpolator& interpolator2, 00862 int num_points1, 00863 int num_points2, 00864 int dim, 00865 const double* param1_start, 00866 const double* param2_start, 00867 const double* data_start); 00868 00880 void gridEvaluator(int num_u, int num_v, 00881 std::vector<double>& points, 00882 std::vector<double>& normals, 00883 std::vector<double>& param_u, 00884 std::vector<double>& param_v, 00885 bool normalize = true) const; 00886 00896 void gridEvaluator(int num_u, int num_v, 00897 std::vector<double>& points, 00898 std::vector<double>& param_u, 00899 std::vector<double>& param_v) const 00900 { 00901 gridEvaluator(num_u, num_v, 00902 points, param_u, param_v, 00903 startparam_u(), 00904 endparam_u(), 00905 startparam_v(), 00906 endparam_v()); 00907 } 00908 00922 void gridEvaluator(int num_u, int num_v, 00923 std::vector<double>& points, 00924 std::vector<double>& param_u, 00925 std::vector<double>& param_v, 00926 double start_u, 00927 double end_u, 00928 double start_v, 00929 double end_v) const; 00930 00938 void gridEvaluator(std::vector<double>& points, 00939 const std::vector<double>& param_u, 00940 const std::vector<double>& param_v) const; 00941 00952 void computeBasis(double param[], 00953 std::vector< double > &basisValues, 00954 std::vector< double > &basisDerivs_u, 00955 std::vector< double > &basisDerivs_v) const; 00956 00973 typedef std::vector<double> Dvector; // for convenience to shorten 00974 typedef std::vector<Dvector> Dmatrix; // function argument lists... 00975 00976 void computeBasisGrid(const Dvector& param_u, 00977 const Dvector& param_v, 00978 Dmatrix& basisValues) const; 00979 00980 00981 void computeBasisGrid(const Dvector& param_u, 00982 const Dvector& param_v, 00983 std::vector<BasisPtsSf>& result) const; 00984 01003 void computeBasisGrid(const Dvector& param_u, 01004 const Dvector& param_v, 01005 Dmatrix& basisValues, 01006 Dmatrix& basisDerivs_u, 01007 Dmatrix& basisDerivs_v) const; 01008 01009 01010 void computeBasisGrid(const Dvector& param_u, 01011 const Dvector& param_v, 01012 std::vector<BasisDerivsSf>& result) const; 01013 01014 01015 void computeBasisGrid(const Dvector& param_u, 01016 const Dvector& param_v, 01017 std::vector<BasisDerivsSf2>& result) const; 01018 01019 01020 // inherited from ParamSurface 01021 virtual double nextSegmentVal(int dir, double par, bool forward, double tol) const; 01022 01023 // Replace one boundary curve of this surface 01024 // Update spline spaces to enable replacement 01025 // NB! Requires both the surface and the new boundary curve to be 01026 // either rational or non-rational 01027 // bd_nmb = 0: umin 01028 // = 1: umax 01029 // = 2: vmin 01030 // = 3: vmax 01031 // Return value: true if a replacement is performe 01032 bool replaceBoundaryCurve(int bd_nmb, boost::shared_ptr<SplineCurve> bd_crv, 01033 bool unify=true); 01034 01036 virtual bool isSpline() const 01037 { 01038 return true; 01039 } 01040 01049 void add(const SplineSurface* other, double tol = 1.0e-10); 01050 01054 void getBoundaryCoefficients(int bd_nmb, std::vector<int>& enumeration) const; 01055 01057 void representAsRational(); 01058 01059 private: 01060 01061 struct degenerate_info 01062 { 01063 bool is_set_; 01064 double tol_; 01065 bool b_, r_, t_, l_; 01066 01067 degenerate_info() 01068 { is_set_ = false; } 01069 }; 01070 01071 // Canonical data 01072 int dim_; 01073 bool rational_; 01074 BsplineBasis basis_u_; 01075 BsplineBasis basis_v_; 01076 std::vector<double> coefs_; 01077 std::vector<double> rcoefs_; 01078 01079 // Generated data 01080 mutable RectDomain domain_; 01081 mutable CurveLoop spatial_boundary_; 01082 mutable degenerate_info degen_; 01083 01084 01085 // Helper functions 01086 void updateCoefsFromRcoefs(); 01087 std::vector<double>& activeCoefs() { return rational_ ? rcoefs_ : coefs_; } 01088 bool normal_not_failsafe(Point& n, double upar, double vpar) const; 01089 bool search_for_normal(bool interval_in_u, 01090 double fixed_parameter, 01091 double interval_start, // normal is not defined here 01092 double interval_end, 01093 Point& normal) const; 01094 01095 // Members new to this class 01096 public: 01097 void pointsGrid(int m1, int m2, int derivs, 01098 const double* basisvals1, 01099 const double* basisvals2, 01100 const int* knotint1, 01101 const int* knotint2, 01102 double* result, 01103 double* normals = 0) const; 01104 private: 01105 01106 // Rewritten pointsGrid, to avoid reformatting results. 01107 // Does not return the derivatives, works only for a 3D non-rational 01108 // spline. 01109 void pointsGridNoDerivs(int m1, int m2, 01110 const double* basisvals1, 01111 const double* basisvals2, 01112 const int* knotint1, 01113 const int* knotint2, 01114 double* result, 01115 double* normals = 0, 01116 bool normalize = true) const; 01117 01118 void accumulateBasis(double* basisvals_u, int uorder, 01119 double* basisvals_v, int vorder, 01120 double* weights, 01121 double* basisValues) const; 01122 01123 void accumulateBasis(double* basisvals_u, int uorder, 01124 double* basisvals_v, int vorder, 01125 double* weights, 01126 double* basisValues, 01127 double* basisDerivs_u, 01128 double* basisDerivs_v) const; 01129 01130 void accumulateBasis(double* basisvals_u, int uorder, 01131 double* basisvals_v, int vorder, 01132 double* weights, 01133 double* basisValues, 01134 double* basisDerivs_u, 01135 double* basisDerivs_v, 01136 double* basisDerivs_uu, 01137 double* basisDerivs_uv, 01138 double* basisDerivs_vv) const; 01139 01140 // Actually computes the closest point. Only difference is the explicit 01141 // robust_seedfind-parameter, which is always true in the virtual 01142 // closestPoint() function, when that function calls closestPointImpl(). 01143 void closestPointImpl(const Point& pt, 01144 double& clo_u, 01145 double& clo_v, 01146 Point& clo_pt, 01147 double& clo_dist, 01148 double epsilon, 01149 const RectDomain* domain_of_interest = NULL, 01150 bool robust_seedfind = true, 01151 double *seed = 0) const; 01152 01153 // Compute the area of the subsurface determined by a subdomain of the parametric space 01154 // Uses cross product of diagonal vectors as an estimation of the area. 01155 double areaDiagonalCross(double tol, 01156 double start_u, double end_u, 01157 double start_v, double end_v) const; 01158 01159 void s1773(const double ppoint[],double aepsge, double estart[],double eend[],double enext[], 01160 double gpos[],int *jstat) const; 01161 01162 void s1773_s9corr(double gd[],double acoef1,double acoef2, 01163 double astart1,double aend1,double astart2,double aend2) const; 01164 01165 void s1773_s9dir(double *cdist,double *cdiff1,double *cdiff2, 01166 double PS[],const double *eval1,std::vector<Point> eval2, 01167 double aepsge, int idim,int *jstat) const; 01168 01169 01170 }; 01171 01172 01173 } // namespace Go 01174 01175 01176 01177 01178 #endif // _SPLINESURFACE_H 01179