![]() |
G+Smo
25.01.0
Geometry + Simulation Modules
|
This class applies hierarchical fitting of parametrized point clouds.
| T | coefficient type |
Inheritance diagram for gsHFitting< d, T >:
Collaboration diagram for gsHFitting< d, T >:Public Member Functions | |
| void | applySmoothing (T lambda, gsSparseMatrix< T > &A_mat) |
| Adds to the matrix A_mat terms for minimization of second derivative, weighted with parameter lambda. | |
| void | assembleSystem (const gsMatrix< T > &points_int, const gsMatrix< T > ¶ms_int, const gsMatrix< T > &points_bdy, const gsMatrix< T > ¶ms_bdy, const index_t &num_basis, const gsSparseMatrix< T > &NNT, gsSparseMatrix< T > &A_tilde, gsMatrix< T > &rhs) |
| assembleSystem: assembles the linear system for the Hybrid Distance Minimization method | |
| void | assembleSystem (gsSparseMatrix< T > &A_mat, gsMatrix< T > &B) |
| Assembles system for the least square fit. | |
| void | blending_weights (const gsSparseMatrix< T > &N_int, const index_t &num_int, const T &mu, const T &sigma, const gsMatrix< T > ¶ms_int, tdm_method method, gsSparseMatrix< T > &NNT) |
| blending_weights: computes the blending weights mu and sigma for the balance mu * PDM + sigma * TDM in the HDM method | |
| void | compute (T lambda=0) |
| compute: Computes the coefficients of the spline geometry via penalized least squares | |
| void | compute_normals (const index_t &num_int, const gsMatrix< T > ¶ms_int, gsSparseMatrix< T > &N_int) |
| compute_normals: Computes the normals of the fitted geometry at the input parameter values params_int | |
| void | compute_tdm (T lambda, T mu, T sigma, const std::vector< index_t > &interpIdx, tdm_method method=hybrid_curvature_pdm_tdm_boundary_pdm) |
| compute_hdm: computes the coefficients of the spline geometry via Hybrid Distance Minimization (HDM) | |
| void | computeApproxError (T &error, int type=0) const |
| void | computeErrors () |
| std::vector< T > | computeErrors (const gsMatrix<> ¶m_values, const gsMatrix<> &points) |
| void | computeMaxNormErrors () |
| gsMatrix< T > | fill_pointWiseErrors (const index_t &num_int, T &max_err_int) |
| vector of length num_int containing all the point-wise errors; store also the max err value in max_err_int | |
| void | get_Error (std::vector< T > &errors, int type=0) const |
| const std::vector< unsigned > & | get_extension () const |
| Returns the chosen cell extension. | |
| const gsBasis< T > & | getBasis () const |
| Returns the basis of the approximation. | |
| std::vector< index_t > | getBoxes (const std::vector< T > &errors, const T threshold) |
| Returns boxes which define refinment area. | |
| T | getRefPercentage () const |
| Return the refinement percentage. | |
| gsHFitting () | |
| Default constructor. | |
| gsHFitting (gsMatrix< T > const ¶m_values, gsMatrix< T > const &points, gsHTensorBasis< d, T > &basis, T refin, const std::vector< unsigned > &extension, T lambda=0) | |
| gsHFitting: Main constructor of the h-fitting class | |
| void | initializeGeometry (const gsMatrix< T > &coefficients, const gsMatrix< T > ¶meters) |
| initializeGeometry: Initializes the fitted geometry with given coefficients and parameters | |
| void | initParametricDomain () |
| Initialize the parametric domain of the point cloud. | |
| gsMatrix< T > | inverse_principal_curvatures (const index_t &num_int, const gsMatrix< T > ¶ms_int) |
| vector of length num_int containing rho = 1/max(c1, c2), where c1, c2 are the principal curvature values computed at every parametric point params_int | |
| bool | is_corner (gsMatrix< T > ¶metric_domain, gsVector< T > ¶meter) |
| check if the given parameter parameter is a corner of the domain parametric_domain | |
| bool | is_point_inside_support (const gsMatrix< T > ¶meter, const gsMatrix< T > &support) |
| bool | is_point_inside_support (const T x, const T y, const gsMatrix< T > &support) |
| check if the given parameter x, y is inside the support support; same as is_point_inside_support, but different input format | |
| bool | is_point_within_cell (const gsMatrix< T > ¶meter, const gsMatrix< T > &element) |
| check if the given parameter parameter is within the cell element | |
| bool | is_point_within_cell (const T x, const T y, const gsMatrix< T > &element) |
| check if the given parameter x, y is within the cell element; same as is_point_within_cell, but different input format | |
| void | iterativeCompute (T const &tolerance, unsigned const &num_iters=10) |
| Computes the least squares fit for a gsBasis. | |
| void | iterativeRefine (int iterations, T tolerance, T err_threshold=-1) |
| iterativeRefine: iteratively refine the basis | |
| T | lambda () const |
| Returns the smoothing weight used in the last fitting. | |
| T | maxPointError () const |
| Returns the maximum point-wise error from the pount cloud (or zero if not fitted) | |
| T | minPointError () const |
| Returns the minimum point-wise error from the pount cloud (or zero if not fitted) | |
| const gsMappedSpline< 2, T > & | mresult () const |
| gives back the computed approximation for multipatch geometry | |
| bool | nextIteration (T tolerance, T err_threshold, const std::vector< boxSide > &fixedSides, index_t maxPcIter=0) |
| Like nextIteration without fixedSides but keeping the values on these sides unchanged throughout the fit. | |
| bool | nextIteration (T tolerance, T err_threshold, index_t maxPcIter=0) |
| nextIteration: perform one iterazion of adaptive refinement with PDM fitting without boundary constraints for iterativeRefine; | |
| bool | nextIteration_pdm (T tolerance, T err_threshold, const std::vector< boxSide > &fixedSides, index_t maxPcIter, const std::vector< index_t > &interpIdx, bool admissibleRef) |
| Like nextIteration_pdm without fixedSides but keeping the values on these sides unchanged throughout the fit. | |
| bool | nextIteration_pdm (T tolerance, T err_threshold, index_t maxPcIter, const std::vector< index_t > &interpIdx, bool admissibleRef=false) |
| nextIteration_pdm: perform one iterazion of adaptive refinement with PDM fitting and with boundary constraints; | |
| bool | nextIteration_tdm (T tolerance, T err_threshold, const std::vector< boxSide > &fixedSides, index_t maxPcIter, T mu, T sigma, const std::vector< index_t > &interpIdx, tdm_method method, bool admissibleRef) |
| Like nextIteration_tdm without fixedSides but keeping the values on these sides unchanged throughout the fit. | |
| bool | nextIteration_tdm (T tolerance, T err_threshold, index_t maxPcIter, T mu, T sigma, const std::vector< index_t > &interpIdx, tdm_method method, bool admissibleRef=false) |
| nextIteration_tdm: perform one iterazion of adaptive refinement with HDM fitting and with boundary constraints; | |
| bool | nextRefinement (T tolerance, T err_threshold, const std::vector< boxSide > &fixedSides, index_t maxPcIter=0, bool admissibleRef=false) |
| Like nextRefinement without fixedSides but keeping the values on these sides unchanged throughout the fit. | |
| bool | nextRefinement (T tolerance, T err_threshold, index_t maxPcIter=0, bool admissibleRef=false) |
| nextRefinement: One step of the refinement of iterativeRefine; | |
| size_t | numPointsBelow (T threshold) const |
| Computes the number of points below the error threshold (or zero if not fitted) | |
| void | parameterCorrection (T accuracy=1e-8, index_t maxIter=10, T tolOrth=1e-6) |
| parameterCorrection: globally apply maxIter steps of parameter correction to the least squares fitted geometry | |
| void | parameterCorrectionSepBoundary_pdm (T accuracy, index_t maxIter, const std::vector< index_t > &sepIndex) |
| parameterCorrectionSepBoundary_pdm: apply maxIter steps of parameter correction for PDM method, separating interior and boundary points | |
| void | parameterCorrectionSepBoundary_tdm (T accuracy, index_t maxIter, T mu, T sigma, const std::vector< index_t > &sepIndex, tdm_method method=hybrid_curvature_pdm_tdm_boundary_pdm) |
| parameterCorrectionSepBoundary_tdm: apply maxIter steps of parameter correction for HDM method, separating interior and boundary points | |
| void | parameterProjectionSepBoundary (T accuracy, const std::vector< index_t > &interpIdx) |
| parameterProjectionSepBoundary: project the points onto the fitted geometry, separating interior and boundary points | |
| const std::vector< T > & | pointWiseErrors () const |
| Return the errors for each point. | |
| gsMatrix< T > | pointWiseErrors (const gsMatrix<> ¶meters, const gsMatrix<> &points) |
| gsMatrix< T > | principal_curvatures (const gsMatrix< T > ¶ms) |
| compute the principal curvatures (c1, c2) at the given parameters params | |
| gsGeometry< T > * | result () const |
| gives back the computed approximation | |
| gsMatrix< T > & | returnParamValues () |
| returns the parameter values | |
| gsMatrix< T > | returnPoints () const |
| returns the points | |
| void | setConstraints (const gsSparseMatrix< T > &lhs, const gsMatrix< T > &rhs) |
| void | setConstraints (const std::vector< boxSide > &fixedSides) |
| void | setConstraints (const std::vector< boxSide > &fixedSides, const std::vector< gsBSpline< T > > &fixedCurves) |
| void | setConstraints (const std::vector< index_t > &indices, const std::vector< gsMatrix< T > > &coefs) |
| void | setExtension (std::vector< unsigned > const &extension) |
| Sets the cell extension. | |
| void | setRefPercentage (double refPercent) |
| Sets the refinement percentage. | |
| void | updateGeometry (gsMatrix< T > coefficients, gsMatrix< T > parameters) |
| updateGeometry: Updates the fitted geometry with new coefficients and parameters | |
Protected Member Functions | |
| virtual void | appendBox (std::vector< index_t > &boxes, std::vector< index_t > &cells, const gsVector< T > ¶meter) |
| Appends a box around parameter to the boxes only if the box is not already in boxes. | |
| void | assembleBlockB (const gsMatrix< T > &points, const gsMatrix< T > ¶ms, index_t num_basis, gsSparseMatrix< T > &result) const |
| Assembles 3xblock collocation matrix. | |
| void | assembleBlockX (const gsMatrix< T > &points, gsMatrix< T > &result) const |
| Assembles the right hand side vectors for PDM/TDM. | |
| gsHBoxContainer< d > | getMarkedHBoxesFromBasis_max (const gsHTensorBasis< d, T > &basis, const std::vector< T > &errors, const gsMatrix< T > ¶meters, T threshold, T extension) |
| getMarkedHBoxesFromBasis_max: returns the markd cells to refine admissibiliy. | |
| T | setRefineThreshold (const std::vector< T > &errors) |
| Automatic set the refinement threshold. | |
Static Protected Member Functions | |
| static void | append (std::vector< index_t > &boxes, const gsVector< index_t > &box) |
| Appends a box to the end of boxes (This function also works for cells) | |
| static bool | isCellAlreadyInserted (const gsVector< index_t, d > &a_cell, const std::vector< index_t > &cells) |
| Checks if a_cell is already inserted in container cells. | |
Protected Attributes | |
| gsFunctionSet< T > * | m_basis |
| Pointer keeping the basis. | |
| gsSparseMatrix< T > | m_constraintsLHS |
| gsMatrix< T > | m_constraintsRHS |
| std::vector< unsigned > | m_ext |
| Size of the extension. | |
| T | m_lambda |
| Smoothing parameter. | |
| T | m_max_error |
| Maximum point-wise error. | |
| T | m_min_error |
| Minimum point-wise error. | |
| gsMappedSpline< 2, T > | m_mresult |
| Pointer keeping the resulting multipatch geometry. | |
| gsMatrix< T > | m_param_values |
| the parameter values of the point cloud | |
| gsMatrix< T > | m_points |
| the points of the point cloud | |
| T | m_ref |
| How many % to refine - 0-1 interval. | |
| gsGeometry< T > * | m_result |
| Pointer keeping the resulting geometry. | |
Private Member Functions | |
| void | extendSystem (gsSparseMatrix< T > &A_mat, gsMatrix< T > &m_B) |
| Extends the system of equations by taking constraints into account. | |
|
inline |
gsHFitting: Main constructor of the h-fitting class
| param_values | a matrix containing the parameter values that parametrize the points |
| points | matrix containing the points to be fitted |
| basis | hiearchical basis to use for fitting |
| refin | percentage of errors to refine (if this strategy is chosen) |
| extension | extension to apply to marked cells |
| lambda | smoothing weight |
|
inherited |
assembleSystem: assembles the linear system for the Hybrid Distance Minimization method
| points_int | interior points |
| params_int | interior parameters |
| points_bdy | boundary points |
| params_bdy | boundary parameters |
| num_basis | dimension of the basis |
| NNT | matrix containing the normals and the blending weights |
| A_tilde | output system matrix |
| rhs | output right-hand side vector |
|
inherited |
blending_weights: computes the blending weights mu and sigma for the balance mu * PDM + sigma * TDM in the HDM method
| N_int | matrix containing the nomals at the parameters params_int |
| num_int | indeces of the interior parameters |
| mu | weight for PDM |
| sigma | weight for TDM |
| params_int | input parameter values |
| method | method for computing the blending weights: constant, based on point-wise error, based on curvature |
| NNT | output matrix containing the normals and the blending weights |
|
inherited |
compute: Computes the coefficients of the spline geometry via penalized least squares
| lambda | smoothing weight |
|
inherited |
compute_normals: Computes the normals of the fitted geometry at the input parameter values params_int
| num_int | index of the input parameter values |
| params_int | input parameter values |
| N_int | matrix containing the normals |
|
inherited |
compute_hdm: computes the coefficients of the spline geometry via Hybrid Distance Minimization (HDM)
| lambda | smoothing weight |
| mu | weight for the PDM |
| sigma | weight for the TDM |
| interpIdx | vector containing the number of interior points and the indices of the boundary points |
| method | method for computing the blending weights |
|
inherited |
Computes the approximation error error of the fitted geometry to the original point cloud type = 0: sum of squares, type = 1: sum of absolute values (l_1 norm)
|
inherited |
Computes the point-wise errors in euclidean norm as well as the max and min errors, and updates the member variables m_pointErrors, m_max_error, m_min_error; different from computeMaxNormErrors(), where the error is computed in inifinity/maximum norm
|
inherited |
Computes min, max and mse errors in euclidean norms between the fitted geometry at the parameters param_values and the input point cloud points it does not update the member variables
|
inherited |
Computes the point-wise errors in infinity/maximum norm as well as the max and min errors, and updates the member variables m_pointErrors, m_max_error, m_min_error; different from computeErrors(), where the error is computed in euclidean norm
|
inherited |
Compute the point-wise errors for each point type = 0: point-wise infinity/maximum norm
|
protected |
getMarkedHBoxesFromBasis_max: returns the markd cells to refine admissibiliy.
| basis | the hierarchical basis from which we extract the elements of the domain |
| error | the pointwise parameter error |
| parameters | the sites on which the point-wise error is computed |
| threshold | the threshold to mark for refinement. |
|
inherited |
initializeGeometry: Initializes the fitted geometry with given coefficients and parameters
| coefficients | the input coefficients |
| parameters | the input parameters |
|
inherited |
check if the given parameter parameter is inside the support support difference with is_point_inside_cell in the inclusion of the left and right interval extremes.
| void iterativeRefine | ( | int | iterations, |
| T | tolerance, | ||
| T | err_threshold = -1 |
||
| ) |
iterativeRefine: iteratively refine the basis
| iterations | maximum number of iterations |
| tolerance | (>=0) if the max error is below the tolerance the refinement stops |
| err_threshold | if non negative all cells with errors bigger than the threshold are refined / If it is equal to -1 the m_ref percentage is used 0 = global refinement |
| bool nextIteration | ( | T | tolerance, |
| T | err_threshold, | ||
| index_t | maxPcIter = 0 |
||
| ) |
nextIteration: perform one iterazion of adaptive refinement with PDM fitting without boundary constraints for iterativeRefine;
| tolerance | (>=0) if the maximum error is below the tolerance the refinement stops; |
| err_threshold | the same as in iterativeRefine. |
| bool nextIteration_pdm | ( | T | tolerance, |
| T | err_threshold, | ||
| index_t | maxPcIter, | ||
| const std::vector< index_t > & | interpIdx, | ||
| bool | admissibleRef = false |
||
| ) |
nextIteration_pdm: perform one iterazion of adaptive refinement with PDM fitting and with boundary constraints;
| tolerance | (>=0) if the maximum error is below the tolerance the refinement stops; |
| err_threshold | the same as in iterativeRefine. |
| interpIdx | is the index of the boundary points to compute with PDM; |
| admissibleRef | if true, the refinement is admissible. |
| bool nextIteration_tdm | ( | T | tolerance, |
| T | err_threshold, | ||
| index_t | maxPcIter, | ||
| T | mu, | ||
| T | sigma, | ||
| const std::vector< index_t > & | interpIdx, | ||
| tdm_method | method, | ||
| bool | admissibleRef = false |
||
| ) |
nextIteration_tdm: perform one iterazion of adaptive refinement with HDM fitting and with boundary constraints;
| tolerance | (>=0) if the maximum error is below the tolerance the refinement stops; |
| err_threshold | the same as in iterative_refine. |
| interpIdx | is the index of the boundary points to compute with PDM; |
| admissibleRef | if true, the refinement is admissible. |
| bool nextRefinement | ( | T | tolerance, |
| T | err_threshold, | ||
| index_t | maxPcIter = 0, |
||
| bool | admissibleRef = false |
||
| ) |
nextRefinement: One step of the refinement of iterativeRefine;
| tolerance | (>=0) if the maximum error is below the tolerance the refinement stops; |
| err_threshold | same as in iterativeRefine; |
| maxPcIter | number of parameter correction steps; |
| admissibleRef | if true, the marking for refinement is admissible. |
|
inherited |
parameterCorrection: globally apply maxIter steps of parameter correction to the least squares fitted geometry
| accuracy | accuracy of the closest point computation |
| maxIter | maximum number of parameter correction steps |
| tolOrth | orthogonality tolerance |
|
inherited |
parameterCorrectionSepBoundary_pdm: apply maxIter steps of parameter correction for PDM method, separating interior and boundary points
| accuracy | accuracy of the closest point computation |
| maxIter | maximum number of parameter correction steps |
| sepIndex | vector containing the number of interior points and the indices of the boundary points |
|
inherited |
parameterCorrectionSepBoundary_tdm: apply maxIter steps of parameter correction for HDM method, separating interior and boundary points
| accuracy | accuracy of the closest point computation |
| maxIter | maximum number of parameter correction steps |
| mu | weight for PDM |
| sigma | weight for TDM |
| sepIndex | vector containing the number of interior points and the indices of the boundary points |
| method | method for computing the blending weights |
|
inherited |
parameterProjectionSepBoundary: project the points onto the fitted geometry, separating interior and boundary points
| accuracy | accuracy of the closest point computation, for the foot-point projection |
| interpIdx | vector containing the number of interior points and the indices of the boundary points |
|
inherited |
Compute point-wise error in euclidean norm between the fitted geometry at the parameters parameters and the input point cloud points similar to computeErrors(), but different input and output format; it does not update the member variables
|
inlineinherited |
Sets constraints that the coefficients of the resulting geometry have to conform to. More precisely, denoting the coefficient vector by x, it enforces lhs * x = rhs.
|
inherited |
Sets constraints in such a way that the previous values at fixedSides of the geometry remain intact.
|
inherited |
Set constraints in such a way that the resulting geometry on each of fixedSides will coincide with the corresponding curve in fixedCurves.
|
inherited |
Sets constraints on that the coefficients of the resulting geometry have to conform to.
| indices | indices (in the coefficient vector) of the prescribed coefficients. |
| coefs | prescribed coefficients. |
|
inherited |
updateGeometry: Updates the fitted geometry with new coefficients and parameters
| coefficients | the new coefficients |
| parameters | the new parameters |
|
protectedinherited |
Left hand-side of the constraints that the coefficients of the resulting geometry have to conform to. This corresponds to matrix D in Prautzch, Boehm, Paluszny: Bezier and B-spline techniques, Section 4.7.
|
protectedinherited |
Right hand-side of the constraints that the coefficients of the resulting geometry have to conform to. This corresponds to vector q in Prautzch, Boehm, Paluszny: Bezier and B-spline techniques, Section 4.7.