1 /* Art_Speaker.cpp
2 *
3 * Copyright (C) 1992-2009,2011,2012,2014-2018,2021 Paul Boersma
4 *
5 * This code is free software; you can redistribute it and/or modify
6 * it under the terms of the GNU General Public License as published by
7 * the Free Software Foundation; either version 2 of the License, or (at
8 * your option) any later version.
9 *
10 * This code is distributed in the hope that it will be useful, but
11 * WITHOUT ANY WARRANTY; without even the implied warranty of
12 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
13 * See the GNU General Public License for more details.
14 *
15 * You should have received a copy of the GNU General Public License
16 * along with this work. If not, see <http://www.gnu.org/licenses/>.
17 */
18
19 #include "Art_Speaker.h"
20
21 #define DLIP 5e-3
22
Art_Speaker_toVocalTract(Art _art,Speaker speaker,double intX[],double intY[],double extX[],double extY[],double * out_bodyX,double * out_bodyY)23 void Art_Speaker_toVocalTract (Art _art, Speaker speaker,
24 double intX [], double intY [], double extX [], double extY [],
25 double *out_bodyX, double *out_bodyY)
26 {
27 double *art = _art -> art;
28 double f = speaker -> relativeSize * 1e-3;
29 struct { double x, y, da; } jaw;
30 struct { double dx, dy; } hyoid;
31 struct { double x, y, r, radius; } body;
32 struct { double x, y, r, a; } teeth;
33 struct { double a; } blade;
34 struct { double x, y, a; } tip;
35 struct { double dx, dy; } lowerLip, upperLip;
36 double HBody_x, HBody_y, HC, Sp, p, a, b;
37
38 /* Determine the position of the hyoid bone (Mermelstein's H). */
39 /* The rest position is a characteristic of the speaker. */
40 /* The stylohyoid muscle pulls the hyoid bone up. */
41 /* The sternohyoid muscle pulls the hyoid bone down. */
42 /* The sphincter muscle pulls the hyoid bone backwards. */
43
44 hyoid.dx = -5 * f * art [(int) kArt_muscle::SPHINCTER];
45 hyoid.dy = 20 * f * (art [(int) kArt_muscle::STYLOHYOID]
46 - art [(int) kArt_muscle::STERNOHYOID]);
47
48 /* The larynx moves up and down with the hyoid bone. */
49 /* Only the lowest point (Mermelstein's K) */
50 /* does not follow completely the horizontal movements. */
51
52 /* Anterior larynx. */
53 intX [1] = -14.0 * f + 0.5 * hyoid.dx; intY [1] = -53.0 * f + hyoid.dy;
54 /* Top of larynx. */
55 intX [2] = -20.0 * f + hyoid.dx; intY [2] = -33.0 * f + hyoid.dy;
56 /* Epiglottis. */
57 intX [3] = -20.0 * f + hyoid.dx; intY [3] = -26.0 * f + hyoid.dy;
58 /* Hyoid bone. */
59 intX [4] = -16.0 * f + hyoid.dx; intY [4] = -26.0 * f + hyoid.dy;
60 /* Posterior larynx. */
61 extX [1] = -22.0 * f + hyoid.dx; extY [1] = -53.0 * f + hyoid.dy;
62 /* Esophagus. */
63 extX [2] = -26.0 * f + hyoid.dx; extY [2] = -40.0 * f + hyoid.dy;
64
65 /* The lower pharynx moves up and down with the hyoid bone. */
66 /* The lower constrictor muscle pulls the rear pharyngeal wall forwards. */
67
68 extX [3] = -34.0 * f + art [(int) kArt_muscle::SPHINCTER] * (5.0 * f);
69 extY [3] = extY [2];
70
71 /* The upper pharynx is fixed at the height of the velum. */
72 /* The upper constrictor muscle pulls the rear pharyngeal wall forwards. */
73
74 extX [5] = -34.0 * f + art [(int) kArt_muscle::SPHINCTER] * (5.0 * f);
75 extY [5] = speaker -> velum.y;
76
77 /* The height of the middle pharynx is in between the lower and upper pharynx. */
78 /* The middle constrictor muscle pulls the rear pharyngeal wall forwards. */
79
80 extX [4] = -34.0 * f + art [(int) kArt_muscle::SPHINCTER] * (5.0 * f);
81 extY [4] = 0.5 * (extY [3] + extY [5]);
82
83 /* Tongue root. */
84
85 jaw.x = -75.0 * f; // position of the condyle
86 jaw.y = 53.0 * f;
87 jaw.da = art [(int) kArt_muscle::MASSETER] * 0.15
88 - art [(int) kArt_muscle::MYLOHYOID] * 0.20;
89 body.x = jaw.x + 81.0 * f * cos (-0.60 + jaw.da)
90 - art [(int) kArt_muscle::STYLOGLOSSUS] * (10.0 * f)
91 + art [(int) kArt_muscle::GENIOGLOSSUS] * (10.0 * f);
92 body.y = jaw.y + 81.0 * f * sin (-0.60 + jaw.da)
93 - art [(int) kArt_muscle::HYOGLOSSUS] * (10.0 * f)
94 + art [(int) kArt_muscle::STYLOGLOSSUS] * (5.0 * f);
95 *out_bodyX = body.x;
96 *out_bodyY = body.y;
97 body.r = hypot (jaw.x - body.x, jaw.y - body.y);
98 body.radius = 20.0 * f;
99 HBody_x = body.x - intX [4];
100 HBody_y = body.y - intY [4];
101 HC = hypot (HBody_x, HBody_y);
102 if (HC <= body.radius) {
103 HC = body.radius;
104 Sp = 0.0; // prevent rounding errors in sqrt (can occur on processors with e.g. 80-bit registers)
105 } else
106 Sp = sqrt (sqr (HC) - sqr (body.radius));
107 a = atan2 (HBody_y, HBody_x);
108 b = asin (body.radius / HC);
109 p = 0.57 * (34.8 * f - Sp);
110 intX [5] = intX [4] + 0.5 * Sp * cos (a + b) - p * sin (a + b);
111 intY [5] = intY [4] + 0.5 * Sp * sin (a + b) + p * cos (a + b);
112 HBody_x = body.x - intX [5];
113 HBody_y = body.y - intY [5];
114 HC = hypot (HBody_x, HBody_y);
115 if (HC <= body.radius) {
116 HC = body.radius;
117 Sp = 0.0;
118 } else
119 Sp = sqrt (sqr (HC) - sqr (body.radius));
120 a = atan2 (HBody_y, HBody_x);
121 b = asin (body.radius / HC);
122 intX [6] = intX [5] + Sp * cos (a + b);
123 intY [6] = intY [5] + Sp * sin (a + b);
124
125 /* Posterior blade. */
126
127 teeth.a = speaker -> lowerTeeth.a + jaw.da;
128 intX [7] = body.x + body.radius * cos (1.73 + teeth.a);
129 intY [7] = body.y + body.radius * sin (1.73 + teeth.a);
130
131 /* Tip. */
132
133 tip.a = (art [(int) kArt_muscle::UPPER_TONGUE]
134 - art [(int) kArt_muscle::LOWER_TONGUE]) * 1.0;
135 blade.a = teeth.a
136 + 0.004 * (body.r - speaker -> neutralBodyDistance) + tip.a;
137 intX [8] = intX [7] + speaker -> tip.length * cos (blade.a);
138 intY [8] = intY [7] + speaker -> tip.length * sin (blade.a);
139
140 /* Jaw. */
141
142 teeth.r = speaker -> lowerTeeth.r;
143 teeth.x = jaw.x + teeth.r * cos (teeth.a);
144 teeth.y = jaw.y + teeth.r * sin (teeth.a);
145 intX [9] = teeth.x + speaker -> teethCavity.dx1;
146 intY [9] = teeth.y + speaker -> teethCavity.dy;
147 intX [10] = teeth.x + speaker -> teethCavity.dx2;
148 intY [10] = intY [9];
149 intX [11] = teeth.x;
150 intY [11] = teeth.y;
151
152 /* Lower lip. */
153
154 lowerLip.dx = speaker -> lowerLip.dx + art [(int) kArt_muscle::ORBICULARIS_ORIS] * 0.02 - 5e-3;
155 lowerLip.dy = speaker -> lowerLip.dy + art [(int) kArt_muscle::ORBICULARIS_ORIS] * 0.01;
156 intX [12] = teeth.x;
157 intY [12] = teeth.y + lowerLip.dy;
158 intX [13] = teeth.x + lowerLip.dx;
159 intY [13] = intY [12];
160
161 /* Velum. */
162
163 extX [6] = speaker -> velum.x;
164 extY [6] = speaker -> velum.y;
165
166 /* Palate. */
167
168 extX [7] = speaker -> alveoli.x;
169 extY [7] = speaker -> alveoli.y;
170 extX [8] = speaker -> upperTeeth.x;
171 extY [8] = speaker -> upperTeeth.y;
172
173 /* Upper lip. */
174
175 upperLip.dx = speaker -> upperLip.dx + art [(int) kArt_muscle::ORBICULARIS_ORIS] * 0.02 - 5e-3;
176 upperLip.dy = speaker -> upperLip.dy - art [(int) kArt_muscle::ORBICULARIS_ORIS] * 0.01;
177 extX [9] = extX [8];
178 extY [9] = extY [8] + upperLip.dy;
179 extX [10] = extX [9] + upperLip.dx;
180 extY [10] = extY [9];
181 extX [11] = extX [10] + 5e-3;
182 extY [11] = extY [10] + DLIP;
183
184 /* Chin. */
185
186 intX [14] = intX [13] + 5e-3;
187 intY [14] = intY [13] - DLIP;
188 intX [15] = intX [11] + 0.5e-2;
189 intY [15] = intY [11] - 3.0e-2;
190 intX [16] = intX [1];
191 intY [16] = intY [1];
192 }
193
Art_Speaker_draw(Art art,Speaker speaker,Graphics g)194 void Art_Speaker_draw (Art art, Speaker speaker, Graphics g) {
195 const double f = speaker -> relativeSize * 1e-3;
196 double intX [1 + 16], intY [1 + 16], extX [1 + 11], extY [1 + 11];
197 double bodyX, bodyY;
198 Art_Speaker_toVocalTract (art, speaker, intX, intY, extX, extY, & bodyX, & bodyY);
199 const Graphics_Viewport previous = Graphics_insetViewport (g, 0.1, 0.9, 0.1, 0.9);
200 Graphics_setWindow (g, -0.05, 0.05, -0.05, 0.05);
201
202 /* Draw inner contour. */
203
204 for (integer i = 1; i <= 5; i ++)
205 Graphics_line (g, intX [i], intY [i], intX [i + 1], intY [i + 1]);
206 Graphics_arc (g, bodyX, bodyY, 20.0 * f,
207 atan2 (intY [7] - bodyY, intX [7] - bodyX) * (180.0 / NUMpi),
208 atan2 (intY [6] - bodyY, intX [6] - bodyX) * (180.0 / NUMpi));
209 for (integer i = 7; i <= 15; i ++)
210 Graphics_line (g, intX [i], intY [i], intX [i + 1], intY [i + 1]);
211
212 /* Draw outer contour. */
213
214 for (integer i = 1; i <= 5; i ++)
215 Graphics_line (g, extX [i], extY [i], extX [i + 1], extY [i + 1]);
216 Graphics_arc (g, 0.0, 0.0, speaker -> palate.radius,
217 speaker -> alveoli.a * (180.0 / NUMpi),
218 speaker -> velum.a * (180.0 / NUMpi));
219 for (integer i = 7; i <= 10; i ++)
220 Graphics_line (g, extX [i], extY [i], extX [i + 1], extY [i + 1]);
221 Graphics_resetViewport (g, previous);
222 }
223
Art_Speaker_fillInnerContour(Art art,Speaker speaker,Graphics g)224 void Art_Speaker_fillInnerContour (Art art, Speaker speaker, Graphics g) {
225 const double f = speaker -> relativeSize * 1e-3;
226 double intX [1 + 16], intY [1 + 16], extX [1 + 11], extY [1 + 11];
227 double x [1 + 16], y [1 + 16];
228 double bodyX, bodyY;
229 Art_Speaker_toVocalTract (art, speaker, intX, intY, extX, extY, & bodyX, & bodyY);
230 const Graphics_Viewport previous = Graphics_insetViewport (g, 0.1, 0.9, 0.1, 0.9);
231 Graphics_setWindow (g, -0.05, 0.05, -0.05, 0.05);
232 for (integer i = 1; i <= 16; i ++) {
233 x [i] = intX [i];
234 y [i] = intY [i];
235 }
236 Graphics_setGrey (g, 0.8);
237 Graphics_fillArea (g, 16, & x [1], & y [1]);
238 Graphics_fillCircle (g, bodyX, bodyY, 20.0 * f);
239 Graphics_setGrey (g, 0.0);
240 Graphics_resetViewport (g, previous);
241 }
242
arcLength(double from,double to)243 static double arcLength (double from, double to) {
244 double result = to - from;
245 while (result > 0.0)
246 result -= 2.0 * NUMpi;
247 while (result < 0.0)
248 result += 2.0 * NUMpi;
249 return result;
250 }
251
252 static int Art_Speaker_meshCount = 27;
253 static double bodyX, bodyY, bodyRadius;
254
toLine(double x,double y,const double intX[],const double intY[],integer i)255 static double toLine (double x, double y, const double intX [], const double intY [], integer i) {
256 integer nearby;
257 if (i == 6) {
258 const double a7 = atan2 (intY [7] - bodyY, intX [7] - bodyX);
259 const double a6 = atan2 (intY [6] - bodyY, intX [6] - bodyX);
260 const double a = atan2 (y - bodyY, x - bodyX);
261 const double da6 = arcLength (a7, a6);
262 const double da = arcLength (a7, a);
263 if (da <= da6)
264 return fabs (hypot (bodyX - x, bodyY - y) - bodyRadius);
265 else
266 nearby = ( arcLength (a7 + 0.5 * da6, a) < NUMpi ? 6 : 7 );
267 } else if ((x - intX [i]) * (intX [i + 1] - intX [i]) +
268 (y - intY [i]) * (intY [i + 1] - intY [i]) < 0) {
269 nearby = i;
270 } else if ((x - intX [i + 1]) * (intX [i] - intX [i + 1]) +
271 (y - intY [i + 1]) * (intY [i] - intY [i + 1]) < 0) {
272 nearby = i + 1;
273 } else {
274 const double boundaryDistance = hypot (intX [i + 1] - intX [i], intY [i + 1] - intY [i]);
275 const double outerProduct = (intX [i] - x) * (intY [i + 1] - intY [i]) - (intY [i] - y) * (intX [i + 1] - intX [i]);
276 return fabs (outerProduct) / boundaryDistance;
277 }
278 return hypot (intX [nearby] - x, intY [nearby] - y);
279 }
280
inside(double x,double y,const double intX[],const double intY[])281 static int inside (double x, double y,
282 const double intX [], const double intY [])
283 {
284 integer up = 0;
285 for (integer i = 1; i <= 16 - 1; i ++)
286 if ((y > intY [i]) != (y > intY [i + 1])) {
287 const double slope = (intX [i + 1] - intX [i]) / (intY [i + 1] - intY [i]);
288 if (x > intX [i] + (y - intY [i]) * slope)
289 up += ( y > intY [i] ? 1 : -1 );
290 }
291 return up != 0 || bodyRadius * bodyRadius > sqr (x - bodyX) + sqr (y - bodyY);
292 }
293
Art_Speaker_meshVocalTract(Art art,Speaker speaker,double xi[],double yi[],double xe[],double ye[],double xmm[],double ymm[],bool closed[])294 void Art_Speaker_meshVocalTract (Art art, Speaker speaker,
295 double xi [], double yi [], double xe [], double ye [],
296 double xmm [], double ymm [], bool closed [])
297 {
298 double f = speaker -> relativeSize * 1e-3;
299 double intX [1 + 16], intY [1 + 16], extX [1 + 11], extY [1 + 11], d_angle;
300 double xm [40], ym [40];
301
302 Art_Speaker_toVocalTract (art, speaker, intX, intY, extX, extY, & bodyX, & bodyY);
303 bodyRadius = 20.0 * f;
304
305 xe [1] = extX [1]; // eq. 5.45
306 ye [1] = extY [1];
307 xe [2] = 0.2 * extX [2] + 0.8 * extX [1];
308 ye [2] = 0.2 * extY [2] + 0.8 * extY [1];
309 xe [3] = 0.6 * extX [2] + 0.4 * extX [1];
310 ye [3] = 0.6 * extY [2] + 0.4 * extY [1];
311 xe [4] = 0.9 * extX [3] + 0.1 * extX [4]; // eq. 5.46
312 ye [4] = 0.9 * extY [3] + 0.1 * extY [4];
313 xe [5] = 0.7 * extX [3] + 0.3 * extX [4];
314 ye [5] = 0.7 * extY [3] + 0.3 * extY [4];
315 xe [6] = 0.5 * extX [3] + 0.5 * extX [4];
316 ye [6] = 0.5 * extY [3] + 0.5 * extY [4];
317 xe [7] = 0.3 * extX [3] + 0.7 * extX [4];
318 ye [7] = 0.3 * extY [3] + 0.7 * extY [4];
319 xe [8] = 0.1 * extX [3] + 0.9 * extX [4];
320 ye [8] = 0.1 * extY [3] + 0.9 * extY [4];
321 xe [9] = 0.9 * extX [4] + 0.1 * extX [5];
322 ye [9] = 0.9 * extY [4] + 0.1 * extY [5];
323 xe [10] = 0.7 * extX [4] + 0.3 * extX [5];
324 ye [10] = 0.7 * extY [4] + 0.3 * extY [5];
325 xe [11] = 0.5 * extX [4] + 0.5 * extX [5];
326 ye [11] = 0.5 * extY [4] + 0.5 * extY [5];
327 xe [12] = 0.3 * extX [4] + 0.7 * extX [5];
328 ye [12] = 0.3 * extY [4] + 0.7 * extY [5];
329 xe [13] = 0.1 * extX [4] + 0.9 * extX [5];
330 ye [13] = 0.1 * extY [4] + 0.9 * extY [5];
331 d_angle = (atan2 (ye [13], xe [13]) - 0.5 * NUMpi) / 6; // eq. 5.47
332 for (integer i = 14; i <= 18; i ++) {
333 double a = 0.5 * NUMpi + (19 - i) * d_angle;
334 xe [i] = speaker -> palate.radius * cos (a);
335 ye [i] = speaker -> palate.radius * sin (a);
336 }
337 xe [19] = 0;
338 ye [19] = speaker -> palate.radius;
339 xe [20] = 0.25 * extX [7];
340 xe [21] = 0.50 * extX [7];
341 xe [22] = 0.75 * extX [7];
342 for (integer i = 20; i <= 22; i ++) {
343 ye [i] = speaker -> palate.radius * sqrt (1.0 - xe [i] * xe [i] /
344 (speaker -> palate.radius * speaker -> palate.radius));
345 }
346 xe [23] = extX [7];
347 ye [23] = extY [7];
348 xe [24] = 0.5 * (extX [7] + extX [8]);
349 ye [24] = 0.5 * (extY [7] + extY [8]);
350 xe [25] = extX [8];
351 ye [25] = extY [8];
352 xe [26] = 0.25 * extX [11] + 0.75 * extX [9];
353 xe [27] = 0.75 * extX [11] + 0.25 * extX [9];
354 ye [26] = extY [10];
355 ye [27] = 0.5 * (extY [10] + extY [11]);
356 for (integer i = 1; i <= 27; i ++) { // every mesh point
357 double minimum = 100000.0;
358 for (integer j = 1; j <= 15 - 1; j ++) { // every internal segment
359 double d = toLine (xe [i], ye [i], intX, intY, j);
360 if (d < minimum) minimum = d;
361 }
362 if (( closed [i] = inside (xe [i], ye [i], intX, intY) ))
363 minimum = - minimum;
364 if (xe [i] >= 0.0) { // vertical line pieces
365 xi [i] = xe [i];
366 yi [i] = ye [i] - minimum;
367 } else if (ye [i] <= 0.0) { // horizontal line pieces
368 xi [i] = xe [i] + minimum;
369 yi [i] = ye [i];
370 } else { // radial line pieces, centre = centre of palate arc
371 double angle = atan2 (ye [i], xe [i]);
372 xi [i] = xe [i] - minimum * cos (angle);
373 yi [i] = ye [i] - minimum * sin (angle);
374 }
375 }
376 for (integer i = 1; i <= Art_Speaker_meshCount; i ++) {
377 xm [i] = 0.5 * (xe [i] + xi [i]);
378 ym [i] = 0.5 * (ye [i] + yi [i]);
379 }
380 for (integer i = 2; i <= Art_Speaker_meshCount; i ++) {
381 xmm [i] = 0.5 * (xm [i - 1] + xm [i]);
382 ymm [i] = 0.5 * (ym [i - 1] + ym [i]);
383 }
384 xmm [1] = 2.0 * xm [1] - xmm [2];
385 ymm [1] = 2.0 * ym [1] - ymm [2];
386 xmm [Art_Speaker_meshCount + 1] = 2.0 * xm [Art_Speaker_meshCount]
387 - xmm [Art_Speaker_meshCount];
388 ymm [Art_Speaker_meshCount + 1] = 2.0 * ym [Art_Speaker_meshCount]
389 - ymm [Art_Speaker_meshCount];
390 }
391
Art_Speaker_drawMesh(Art art,Speaker speaker,Graphics graphics)392 void Art_Speaker_drawMesh (Art art, Speaker speaker, Graphics graphics) {
393 double xi [40], yi [40], xe [40], ye [40], xmm [40], ymm [40];
394 bool closed [40];
395 int oldLineType = Graphics_inqLineType (graphics);
396 Art_Speaker_meshVocalTract (art, speaker, xi, yi, xe, ye, xmm, ymm, closed);
397 Graphics_Viewport previous = Graphics_insetViewport (graphics, 0.1, 0.9, 0.1, 0.9); // must be square
398 Graphics_setWindow (graphics, -0.05, 0.05, -0.05, 0.05);
399
400 /* Mesh lines. */
401 for (integer i = 1; i <= Art_Speaker_meshCount; i ++)
402 Graphics_line (graphics, xi [i], yi [i], xe [i], ye [i]);
403
404 /* Radii. */
405 Graphics_setLineType (graphics, Graphics_DOTTED);
406 for (integer i = 1; i <= Art_Speaker_meshCount; i ++)
407 if (xe [i] <= 0.0 && ye [i] >= 0.0)
408 Graphics_line (graphics, 0.0, 0.0, 0.9 * xi [i], 0.9 * yi [i]);
409 Graphics_setLineType (graphics, oldLineType);
410
411 /* Lengths. */
412 for (integer i = 1; i <= Art_Speaker_meshCount; i ++)
413 Graphics_line (graphics, xmm [i], ymm [i], xmm [i + 1], ymm [i + 1]);
414
415 for (integer i = 1; i <= Art_Speaker_meshCount + 1; i ++)
416 Graphics_speckle (graphics, xmm [i], ymm [i]);
417 Graphics_setTextAlignment (graphics, Graphics_LEFT, Graphics_HALF);
418 Graphics_text (graphics, 0.0, 0.0, U"O"); // origin
419 Graphics_resetViewport (graphics, previous);
420 }
421
422 /* End of file Art_Speaker.cpp */
423