1 /*
2 * Copyright (C) 1998, 2000-2007, 2010, 2011, 2012, 2013 SINTEF ICT,
3 * Applied Mathematics, Norway.
4 *
5 * Contact information: E-mail: tor.dokken@sintef.no
6 * SINTEF ICT, Department of Applied Mathematics,
7 * P.O. Box 124 Blindern,
8 * 0314 Oslo, Norway.
9 *
10 * This file is part of SISL.
11 *
12 * SISL is free software: you can redistribute it and/or modify
13 * it under the terms of the GNU Affero General Public License as
14 * published by the Free Software Foundation, either version 3 of the
15 * License, or (at your option) any later version.
16 *
17 * SISL is distributed in the hope that it will be useful,
18 * but WITHOUT ANY WARRANTY; without even the implied warranty of
19 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
20 * GNU Affero General Public License for more details.
21 *
22 * You should have received a copy of the GNU Affero General Public
23 * License along with SISL. If not, see
24 * <http://www.gnu.org/licenses/>.
25 *
26 * In accordance with Section 7(b) of the GNU Affero General Public
27 * License, a covered work must retain the producer line in every data
28 * file that is created or manipulated using SISL.
29 *
30 * Other Usage
31 * You can be released from the requirements of the license by purchasing
32 * a commercial license. Buying such a license is mandatory as soon as you
33 * develop commercial activities involving the SISL library without
34 * disclosing the source code of your own applications.
35 *
36 * This file may be used in accordance with the terms contained in a
37 * written agreement between you and SINTEF ICT.
38 */
39
40 #include "sisl-copyright.h"
41 #define SH6CVVERT
42
43 #include "sislP.h"
44
45
46 #if defined(SISLNEEDPROTOTYPES)
47 void
sh6cvvert(SISLCurve * pc1,SISLCurve * pc2,double * cpar1,double * cpar2)48 sh6cvvert(SISLCurve *pc1, SISLCurve *pc2, double *cpar1, double *cpar2)
49 #else
50 void sh6cvvert(pc1,pc2,cpar1,cpar2)
51 SISLCurve *pc1;
52 SISLCurve *pc2;
53 double *cpar1;
54 double *cpar2;
55 #endif
56 /*
57 *********************************************************************
58 *
59 * PURPOSE : Estimate the parameter values the closest vertices of
60 * two B-spline curves
61 *
62 *
63 *
64 * INPUT : pc1 - Pointer to first curve.
65 * pc1 - Pointer to second curve.
66 *
67 *
68 * OUTPUT : cpar1 - Parameter value of closest vertex of the 1. curve.
69 * cpar2 - Parameter value of closest vertex of the 2. curve.
70 *
71 *
72 * METHOD : Find the vertices of the two curves that are closest
73 * to each other. Regard these vertices as Schoenberg
74 * points, and compute the parameter values.
75 *
76 *
77 * REFERENCES :
78 *
79 *
80 * USE :
81 *
82 *-
83 * CALLS : s6dist - Distance between two points.
84 *
85 *
86 * WRITTEN BY : Vibeke Skytt, SI, 09.00.
87 *
88 *********************************************************************
89 */
90 {
91 int ki,kj, kh; /* Counters. */
92 int kdim = pc1->idim; /* Dimension of geometry space. */
93 int kminc1; /* Number of closest vertex of 1. curve. */
94 int kminc2; /* Number of closest vertex of 2. curve. */
95 int kn1 = pc1->in; /* Number of coefficients of 1. curve. */
96 int kn2 = pc2->in; /* Number of coefficients of 2. curve. */
97 int kk1 = pc1->ik; /* Order of 1. curve. */
98 int kk2 = pc2->ik; /* Order of 2. curve. */
99 double tdist; /* Distance. */
100 double tmin = HUGE; /* Minimum distance. */
101 double tpar; /* Used to compute parameter values. */
102 double *s1,*s2; /* Pointers into arrays. */
103
104 /* Find position of closest vertices. */
105
106 for (s1=pc1->ecoef, ki=0; ki<kn1; s1+=kdim, ki++)
107 for (s2=pc2->ecoef, kj=0; kj<kn2; s2+=kdim, kj++)
108 {
109 for (tdist=0.0, kh=kdim-1; kh>=0; kh--)
110 tdist += (s2[kh]-s1[kh])*(s2[kh]-s1[kh]);
111 // tdist = s6dist(s1,s2,kdim);
112 if (tdist < tmin)
113 {
114 tmin = tdist;
115 kminc1 = ki;
116 kminc2 = kj;
117 }
118 }
119
120 /* Estimate parameter values of vertices. */
121
122 for (ki=kminc1+1, s1=pc1->et+ki, tpar=0.0;
123 ki<kminc1+kk1; tpar+=(*s1), s1++, ki++);
124 *cpar1 = tpar/(double)(kk1-1);
125
126 for (ki=kminc2+1, s1=pc2->et+ki, tpar=0.0;
127 ki<kminc2+kk2; tpar+=(*s1), s1++, ki++);
128 *cpar2 = tpar/(double)(kk2-1);
129
130 return;
131 }
132