Branch data Line data Source code
1 : : /* -*- Mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */
2 : : /*
3 : : * This file is part of the LibreOffice project.
4 : : *
5 : : * This Source Code Form is subject to the terms of the Mozilla Public
6 : : * License, v. 2.0. If a copy of the MPL was not distributed with this
7 : : * file, You can obtain one at http://mozilla.org/MPL/2.0/.
8 : : *
9 : : * This file incorporates work covered by the following license notice:
10 : : *
11 : : * Licensed to the Apache Software Foundation (ASF) under one or more
12 : : * contributor license agreements. See the NOTICE file distributed
13 : : * with this work for additional information regarding copyright
14 : : * ownership. The ASF licenses this file to you under the Apache
15 : : * License, Version 2.0 (the "License"); you may not use this file
16 : : * except in compliance with the License. You may obtain a copy of
17 : : * the License at http://www.apache.org/licenses/LICENSE-2.0 .
18 : : */
19 : :
20 : : #ifndef _HOMMATRIX_TEMPLATE_HXX
21 : : #define _HOMMATRIX_TEMPLATE_HXX
22 : :
23 : : #include <sal/types.h>
24 : : #include <basegfx/numeric/ftools.hxx>
25 : : #include <math.h>
26 : : #include <string.h>
27 : :
28 : : namespace basegfx
29 : : {
30 : : namespace internal
31 : : {
32 : :
33 : 21455855 : inline double implGetDefaultValue(sal_uInt16 nRow, sal_uInt16 nColumn)
34 : : {
35 [ + + ]: 21455855 : if(nRow == nColumn)
36 : 8329132 : return 1.0;
37 : 21455855 : return 0.0;
38 : : }
39 : :
40 : : template < unsigned int _RowSize > class ImplMatLine
41 : : {
42 : : enum { RowSize = _RowSize };
43 : :
44 : : double mfValue[RowSize];
45 : :
46 : : public:
47 : 4286788 : ImplMatLine()
48 : : {
49 : 4286788 : }
50 : :
51 : 5463 : explicit ImplMatLine(sal_uInt16 nRow, ImplMatLine< RowSize >* pToBeCopied = 0L)
52 : : {
53 [ + + ]: 5463 : if(pToBeCopied)
54 : : {
55 : 2656 : memcpy(&mfValue, pToBeCopied, sizeof(double) * RowSize);
56 : : }
57 : : else
58 : : {
59 [ + + ]: 13970 : for(sal_uInt16 a(0); a < RowSize; a++)
60 : : {
61 : 11163 : mfValue[a] = implGetDefaultValue(nRow, a);
62 : : }
63 : : }
64 : 5463 : }
65 : :
66 : 65909029 : double get(sal_uInt16 nColumn) const
67 : : {
68 : 65909029 : return mfValue[nColumn];
69 : : }
70 : :
71 : 9966253 : void set(sal_uInt16 nColumn, const double& rValue)
72 : : {
73 : 9966253 : mfValue[nColumn] = rValue;
74 : 9966253 : }
75 : : };
76 : :
77 : : template < unsigned int _RowSize > class ImplHomMatrixTemplate
78 : : {
79 : : enum { RowSize = _RowSize };
80 : :
81 : : ImplMatLine< RowSize > maLine[RowSize - 1];
82 : : ImplMatLine< RowSize >* mpLine;
83 : :
84 : : public:
85 : : // Is last line used?
86 : 4678534 : bool isLastLineDefault() const
87 : : {
88 [ + + ]: 4678534 : if(!mpLine)
89 : 4673342 : return true;
90 : :
91 [ + - ]: 6744 : for(sal_uInt16 a(0); a < RowSize; a++)
92 : : {
93 : 6744 : const double fDefault(implGetDefaultValue((RowSize - 1), a));
94 : 6744 : const double fLineValue(mpLine->get(a));
95 : :
96 [ + + ]: 6744 : if(!::basegfx::fTools::equal(fDefault, fLineValue))
97 : : {
98 : 5192 : return false;
99 : : }
100 : : }
101 : :
102 : : // reset last line, it equals default
103 : 0 : delete ((ImplHomMatrixTemplate< RowSize >*)this)->mpLine;
104 : 0 : ((ImplHomMatrixTemplate< RowSize >*)this)->mpLine = 0L;
105 : :
106 : 4678534 : return true;
107 : : }
108 : :
109 : 193808 : ImplHomMatrixTemplate()
110 [ + + ]: 638185 : : mpLine(0L)
111 : : {
112 : : // complete initialization with identity matrix, all lines
113 : : // were initialized with a trailing 1 followed by 0's.
114 [ + + ]: 638185 : for(sal_uInt16 a(0); a < RowSize-1; a++)
115 : : {
116 [ + + ]: 1947791 : for(sal_uInt16 b(0); b < RowSize; b++)
117 : 1503414 : maLine[a].set(b, implGetDefaultValue(a, b) );
118 : : }
119 : 193808 : }
120 : :
121 : 1807631 : ImplHomMatrixTemplate(const ImplHomMatrixTemplate& rToBeCopied)
122 [ + + ]: 5650042 : : mpLine(0L)
123 : : {
124 : : // complete initialization using copy
125 [ + + ]: 5650042 : for(sal_uInt16 a(0); a < (RowSize - 1); a++)
126 : : {
127 : 3842411 : memcpy(&maLine[a], &rToBeCopied.maLine[a], sizeof(ImplMatLine< RowSize >));
128 : : }
129 : :
130 [ + + ]: 1807631 : if(rToBeCopied.mpLine)
131 : : {
132 : 2656 : mpLine = new ImplMatLine< RowSize >((RowSize - 1), rToBeCopied.mpLine);
133 : : }
134 : 1807631 : }
135 : :
136 : 2000957 : ~ImplHomMatrixTemplate()
137 : : {
138 [ + + ]: 2000957 : if(mpLine)
139 : : {
140 : 5463 : delete mpLine;
141 : : }
142 : 2000957 : }
143 : :
144 : 62345 : sal_uInt16 getEdgeLength() const { return RowSize; }
145 : :
146 : 79541234 : double get(sal_uInt16 nRow, sal_uInt16 nColumn) const
147 : : {
148 [ + + ]: 79541234 : if(nRow < (RowSize - 1))
149 : : {
150 : 65807007 : return maLine[nRow].get(nColumn);
151 : : }
152 : :
153 [ + + ]: 13734227 : if(mpLine)
154 : : {
155 : 86756 : return mpLine->get(nColumn);
156 : : }
157 : :
158 : 79541234 : return implGetDefaultValue((RowSize - 1), nColumn);
159 : : }
160 : :
161 : 11147945 : void set(sal_uInt16 nRow, sal_uInt16 nColumn, const double& rValue)
162 : : {
163 [ + + ]: 11147945 : if(nRow < (RowSize - 1))
164 : : {
165 : 8442678 : maLine[nRow].set(nColumn, rValue);
166 : : }
167 [ + + ]: 2705267 : else if(mpLine)
168 : : {
169 : 17354 : mpLine->set(nColumn, rValue);
170 : : }
171 : : else
172 : : {
173 : 2687913 : const double fDefault(implGetDefaultValue((RowSize - 1), nColumn));
174 : :
175 [ + + ]: 2687913 : if(!::basegfx::fTools::equal(fDefault, rValue))
176 : : {
177 [ + - ]: 2807 : mpLine = new ImplMatLine< RowSize >((RowSize - 1), 0L);
178 : 2687913 : mpLine->set(nColumn, rValue);
179 : : }
180 : : }
181 : 11147945 : }
182 : :
183 : 667074 : void testLastLine()
184 : : {
185 [ + + ]: 667074 : if(mpLine)
186 : : {
187 : 3190 : bool bNecessary(false);
188 : :
189 [ + + ][ + - ]: 11712 : for(sal_uInt16 a(0);!bNecessary && a < RowSize; a++)
[ + + ]
190 : : {
191 : 8522 : const double fDefault(implGetDefaultValue((RowSize - 1), a));
192 : 8522 : const double fLineValue(mpLine->get(a));
193 : :
194 [ + + ]: 8522 : if(!::basegfx::fTools::equal(fDefault, fLineValue))
195 : : {
196 : 3190 : bNecessary = true;
197 : : }
198 : : }
199 : :
200 [ - + ]: 3190 : if(!bNecessary)
201 : : {
202 : 0 : delete mpLine;
203 : 0 : mpLine = 0L;
204 : : }
205 : : }
206 : 667074 : }
207 : :
208 : : // Left-upper decompositon
209 : 63622 : bool ludcmp(sal_uInt16 nIndex[], sal_Int16& nParity)
210 : : {
211 : : double fBig, fSum, fDum;
212 : : double fStorage[RowSize];
213 : : sal_uInt16 a, b, c;
214 : :
215 : : // #i30874# Initialize nAMax (compiler warns)
216 : 63622 : sal_uInt16 nAMax = 0;
217 : :
218 : 63622 : nParity = 1;
219 : :
220 : : // Calc the max of each line. If a line is empty,
221 : : // stop immediately since matrix is not invertible then.
222 [ + + ]: 260965 : for(a = 0; a < RowSize; a++)
223 : : {
224 : 197343 : fBig = 0.0;
225 : :
226 [ + + ]: 815280 : for(b = 0; b < RowSize; b++)
227 : : {
228 [ + - ]: 617937 : double fTemp(fabs(get(a, b)));
229 : :
230 [ + + ]: 617937 : if(::basegfx::fTools::more(fTemp, fBig))
231 : : {
232 : 258289 : fBig = fTemp;
233 : : }
234 : : }
235 : :
236 [ - + ]: 197343 : if(::basegfx::fTools::equalZero(fBig))
237 : : {
238 : 0 : return false;
239 : : }
240 : :
241 : 197343 : fStorage[a] = 1.0 / fBig;
242 : : }
243 : :
244 : : // start normalizing
245 [ + + ]: 260965 : for(b = 0; b < RowSize; b++)
246 : : {
247 [ + + ]: 407640 : for(a = 0; a < b; a++)
248 : : {
249 [ + - ]: 210297 : fSum = get(a, b);
250 : :
251 [ + + ]: 293350 : for(c = 0; c < a; c++)
252 : : {
253 [ + - ][ + - ]: 83053 : fSum -= get(a, c) * get(c, b);
254 : : }
255 : :
256 [ + - ]: 210297 : set(a, b, fSum);
257 : : }
258 : :
259 : 197343 : fBig = 0.0;
260 : :
261 [ + + ]: 604983 : for(a = b; a < RowSize; a++)
262 : : {
263 [ + - ]: 407640 : fSum = get(a, b);
264 : :
265 [ + + ]: 700990 : for(c = 0; c < b; c++)
266 : : {
267 [ + - ][ + - ]: 293350 : fSum -= get(a, c) * get(c, b);
268 : : }
269 : :
270 [ + - ]: 407640 : set(a, b, fSum);
271 : 407640 : fDum = fStorage[a] * fabs(fSum);
272 : :
273 [ + + ]: 407640 : if(::basegfx::fTools::moreOrEqual(fDum, fBig))
274 : : {
275 : 198143 : fBig = fDum;
276 : 198143 : nAMax = a;
277 : : }
278 : : }
279 : :
280 [ + + ]: 197343 : if(b != nAMax)
281 : : {
282 [ + + ]: 3948 : for(c = 0; c < RowSize; c++)
283 : : {
284 [ + - ]: 3148 : fDum = get(nAMax, c);
285 [ + - ][ + - ]: 3148 : set(nAMax, c, get(b, c));
286 [ + - ]: 3148 : set(b, c, fDum);
287 : : }
288 : :
289 : 800 : nParity = -nParity;
290 : 800 : fStorage[nAMax] = fStorage[b];
291 : : }
292 : :
293 : 197343 : nIndex[b] = nAMax;
294 : :
295 : : // here the failure of precision occurs
296 [ + - ]: 197343 : const double fValBB(fabs(get(b, b)));
297 : :
298 [ - + ]: 197343 : if(::basegfx::fTools::equalZero(fValBB))
299 : : {
300 : 0 : return false;
301 : : }
302 : :
303 [ + + ]: 197343 : if(b != (RowSize - 1))
304 : : {
305 [ + - ]: 133721 : fDum = 1.0 / get(b, b);
306 : :
307 [ + + ]: 344018 : for(a = b + 1; a < RowSize; a++)
308 : : {
309 [ + - ][ + - ]: 210297 : set(a, b, get(a, b) * fDum);
310 : : }
311 : : }
312 : : }
313 : :
314 : 63622 : return true;
315 : : }
316 : :
317 : 192235 : void lubksb(const sal_uInt16 nIndex[], double fRow[]) const
318 : : {
319 : : sal_uInt16 b, ip;
320 : 192235 : sal_Int16 a, a2 = -1;
321 : : double fSum;
322 : :
323 [ + + ]: 789740 : for(a = 0; a < RowSize; a++)
324 : : {
325 : 597505 : ip = nIndex[a];
326 : 597505 : fSum = fRow[ip];
327 : 597505 : fRow[ip] = fRow[a];
328 : :
329 [ + + ]: 597505 : if(a2 >= 0)
330 : : {
331 [ + + ]: 483215 : for(b = a2; b < a; b++)
332 : : {
333 [ + - ]: 280580 : fSum -= get(a, b) * fRow[b];
334 : : }
335 : : }
336 [ + + ]: 394870 : else if(!::basegfx::fTools::equalZero(fSum))
337 : : {
338 : 192235 : a2 = a;
339 : : }
340 : :
341 : 597505 : fRow[a] = fSum;
342 : : }
343 : :
344 [ + + ]: 789740 : for(a = (RowSize - 1); a >= 0; a--)
345 : : {
346 : 597505 : fSum = fRow[a];
347 : :
348 [ + + ]: 1236610 : for(b = a + 1; b < RowSize; b++)
349 : : {
350 [ + - ]: 639105 : fSum -= get(a, b) * fRow[b];
351 : : }
352 : :
353 [ + - ]: 597505 : const double fValueAA(get(a, a));
354 : :
355 [ + - ]: 597505 : if(!::basegfx::fTools::equalZero(fValueAA))
356 : : {
357 [ + - ]: 597505 : fRow[a] = fSum / get(a, a);
358 : : }
359 : : }
360 : 192235 : }
361 : :
362 : 2391143 : bool isIdentity() const
363 : : {
364 : : // last line needs no testing if not existing
365 : : const sal_uInt16 nMaxLine(
366 [ + + ]: 2391143 : sal::static_int_cast<sal_uInt16>(mpLine ? RowSize : (RowSize - 1)) );
367 : :
368 [ + + ]: 2487212 : for(sal_uInt16 a(0); a < nMaxLine; a++)
369 : : {
370 [ + + ]: 3089192 : for(sal_uInt16 b(0); b < RowSize; b++)
371 : : {
372 : 2993123 : const double fDefault(implGetDefaultValue(a, b));
373 [ + - ]: 2993123 : const double fValueAB(get(a, b));
374 : :
375 [ + + ]: 2993123 : if(!::basegfx::fTools::equal(fDefault, fValueAB))
376 : : {
377 : 2375641 : return false;
378 : : }
379 : : }
380 : : }
381 : :
382 : 2391143 : return true;
383 : : }
384 : :
385 : 0 : bool isInvertible() const
386 : : {
387 [ # # ]: 0 : ImplHomMatrixTemplate aWork(*this);
388 : : sal_uInt16 nIndex[RowSize];
389 : : sal_Int16 nParity;
390 : :
391 [ # # ]: 0 : return aWork.ludcmp(nIndex, nParity);
392 : : }
393 : :
394 : : bool isNormalized() const
395 : : {
396 : : if(!mpLine)
397 : : return true;
398 : :
399 : : const double fHomValue(get((RowSize - 1), (RowSize - 1)));
400 : :
401 : : if(::basegfx::fTools::equalZero(fHomValue))
402 : : {
403 : : return true;
404 : : }
405 : :
406 : : const double fOne(1.0);
407 : :
408 : : if(::basegfx::fTools::equal(fOne, fHomValue))
409 : : {
410 : : return true;
411 : : }
412 : :
413 : : return false;
414 : : }
415 : :
416 : 62345 : void doInvert(const ImplHomMatrixTemplate& rWork, const sal_uInt16 nIndex[])
417 : : {
418 : : double fArray[RowSize];
419 : :
420 [ + + ]: 254580 : for(sal_uInt16 a(0); a < RowSize; a++)
421 : : {
422 : : // prepare line
423 : : sal_uInt16 b;
424 [ + + ]: 789740 : for( b = 0; b < RowSize; b++)
425 : : {
426 : 597505 : fArray[b] = implGetDefaultValue(a, b);
427 : : }
428 : :
429 : : // expand line
430 [ + - ]: 192235 : rWork.lubksb(nIndex, fArray);
431 : :
432 : : // copy line transposed to this matrix
433 [ + + ]: 789740 : for( b = 0; b < RowSize; b++)
434 : : {
435 [ + - ]: 597505 : set(b, a, fArray[b]);
436 : : }
437 : : }
438 : :
439 : : // evtl. get rid of last matrix line
440 : 62345 : testLastLine();
441 : 62345 : }
442 : :
443 : : void doNormalize()
444 : : {
445 : : if(mpLine)
446 : : {
447 : : const double fHomValue(get((RowSize - 1), (RowSize - 1)));
448 : :
449 : : for(sal_uInt16 a(0); a < RowSize; a++)
450 : : {
451 : : for(sal_uInt16 b(0); b < RowSize; b++)
452 : : {
453 : : set(a, b, get(a, b) / fHomValue);
454 : : }
455 : : }
456 : :
457 : : // evtl. get rid of last matrix line
458 : : testLastLine();
459 : : }
460 : : }
461 : :
462 : 1277 : double doDeterminant() const
463 : : {
464 [ + - ]: 1277 : ImplHomMatrixTemplate aWork(*this);
465 : : sal_uInt16 nIndex[RowSize];
466 : : sal_Int16 nParity;
467 : 1277 : double fRetval(0.0);
468 : :
469 [ + - ][ + - ]: 1277 : if(aWork.ludcmp(nIndex, nParity))
470 : : {
471 : 1277 : fRetval = (double)nParity;
472 : :
473 : : // last line needs no multiply if not existing; default value would be 1.
474 : : const sal_uInt16 nMaxLine(
475 [ - + ]: 1277 : sal::static_int_cast<sal_uInt16>(aWork.mpLine ? RowSize : (RowSize - 1)) );
476 : :
477 [ + + ]: 5108 : for(sal_uInt16 a(0); a < nMaxLine; a++)
478 : : {
479 [ + - ]: 3831 : fRetval *= aWork.get(a, a);
480 : : }
481 : : }
482 : :
483 : 1277 : return fRetval;
484 : : }
485 : :
486 : : double doTrace() const
487 : : {
488 : : double fTrace = (mpLine) ? 0.0 : 1.0;
489 : : const sal_uInt16 nMaxLine(
490 : : sal::static_int_cast<sal_uInt16>(mpLine ? RowSize : (RowSize - 1)) );
491 : :
492 : : for(sal_uInt16 a(0); a < nMaxLine; a++)
493 : : {
494 : : fTrace += get(a, a);
495 : : }
496 : :
497 : : return fTrace;
498 : : }
499 : :
500 : : void doTranspose()
501 : : {
502 : : for(sal_uInt16 a(0); a < (RowSize - 1); a++)
503 : : {
504 : : for(sal_uInt16 b(a + 1); b < RowSize; b++)
505 : : {
506 : : const double fTemp(get(a, b));
507 : : set(a, b, get(b, a));
508 : : set(b, a, fTemp);
509 : : }
510 : : }
511 : :
512 : : testLastLine();
513 : : }
514 : :
515 : 0 : void doAddMatrix(const ImplHomMatrixTemplate& rMat)
516 : : {
517 [ # # ]: 0 : for(sal_uInt16 a(0); a < RowSize; a++)
518 : : {
519 [ # # ]: 0 : for(sal_uInt16 b(0); b < RowSize; b++)
520 : : {
521 [ # # ]: 0 : set(a, b, get(a, b) + rMat.get(a, b));
522 : : }
523 : : }
524 : :
525 : 0 : testLastLine();
526 : 0 : }
527 : :
528 : 0 : void doSubMatrix(const ImplHomMatrixTemplate& rMat)
529 : : {
530 [ # # ]: 0 : for(sal_uInt16 a(0); a < RowSize; a++)
531 : : {
532 [ # # ]: 0 : for(sal_uInt16 b(0); b < RowSize; b++)
533 : : {
534 [ # # ]: 0 : set(a, b, get(a, b) - rMat.get(a, b));
535 : : }
536 : : }
537 : :
538 : 0 : testLastLine();
539 : 0 : }
540 : :
541 : 0 : void doMulMatrix(const double& rfValue)
542 : : {
543 [ # # ]: 0 : for(sal_uInt16 a(0); a < RowSize; a++)
544 : : {
545 [ # # ]: 0 : for(sal_uInt16 b(0); b < RowSize; b++)
546 : : {
547 [ # # ]: 0 : set(a, b, get(a, b) * rfValue);
548 : : }
549 : : }
550 : :
551 : 0 : testLastLine();
552 : 0 : }
553 : :
554 : 604729 : void doMulMatrix(const ImplHomMatrixTemplate& rMat)
555 : : {
556 : : // create a copy as source for the original values
557 [ + - ]: 604729 : const ImplHomMatrixTemplate aCopy(*this);
558 : :
559 : : // TODO: maybe optimize cases where last line is [0 0 1].
560 : :
561 : 604729 : double fValue(0.0);
562 : :
563 [ + + ]: 2510103 : for(sal_uInt16 a(0); a < RowSize; ++a)
564 : : {
565 [ + + ]: 7986244 : for(sal_uInt16 b(0); b < RowSize; ++b)
566 : : {
567 : 6080870 : fValue = 0.0;
568 : :
569 [ + + ]: 25782472 : for(sal_uInt16 c(0); c < RowSize; ++c)
570 [ + - ][ + - ]: 19701602 : fValue += aCopy.get(c, b) * rMat.get(a, c);
571 : :
572 [ + - ]: 6080870 : set(a, b, fValue);
573 : : }
574 : : }
575 : :
576 : 604729 : testLastLine();
577 : 604729 : }
578 : :
579 : 63434 : bool isEqual(const ImplHomMatrixTemplate& rMat) const
580 : : {
581 : : const sal_uInt16 nMaxLine(
582 [ + + ][ - + ]: 63434 : sal::static_int_cast<sal_uInt16>((mpLine || rMat.mpLine) ? RowSize : (RowSize - 1)) );
583 : :
584 [ + + ]: 180885 : for(sal_uInt16 a(0); a < nMaxLine; a++)
585 : : {
586 [ + + ]: 490242 : for(sal_uInt16 b(0); b < RowSize; b++)
587 : : {
588 [ + - ]: 372791 : const double fValueA(get(a, b));
589 [ + - ]: 372791 : const double fValueB(rMat.get(a, b));
590 : :
591 [ + + ]: 372791 : if(!::basegfx::fTools::equal(fValueA, fValueB))
592 : : {
593 : 6271 : return false;
594 : : }
595 : : }
596 : : }
597 : :
598 : 63434 : return true;
599 : : }
600 : : };
601 : :
602 : : } // namespace internal
603 : : } // namespace basegfx
604 : :
605 : : #endif /* _HOMMATRIX_TEMPLATE_HXX */
606 : :
607 : : /* vim:set shiftwidth=4 softtabstop=4 expandtab: */
|
