@@ -272,182 +272,4 @@ namespace RadeonProRender
272
272
#undef m
273
273
}
274
274
275
- // /
276
- // / Matrix decomposition based on: http://tog.acm.org/resources/GraphicsGems/gemsii/unmatrix.c
277
- // /
278
- inline int
279
- decompose (matrix const & mat,
280
- float3& translation,
281
- float3& scale,
282
- float3& shear,
283
- float3& rotation,
284
- float4& perspective)
285
- {
286
- int i, j;
287
- matrix locmat;
288
- matrix pmat, invpmat, tinvpmat;
289
- /* Vector4 type and functions need to be added to the common set. */
290
- float4 prhs, psol;
291
- float4 row[3 ], pdum3;
292
-
293
- locmat = mat;
294
- /* Normalize the matrix. */
295
- if (locmat.m33 == 0 )
296
- return 0 ;
297
- for (i = 0 ; i<4 ; i++)
298
- for (j = 0 ; j<4 ; j++)
299
- locmat.m [i][j] /= locmat.m33 ;
300
- /* pmat is used to solve for perspective, but it also provides
301
- * an easy way to test for singularity of the upper 3x3 component.
302
- */
303
-
304
- pmat = locmat;
305
- for (i = 0 ; i<3 ; i++)
306
- pmat.m [i][3 ] = 0 ;
307
- pmat.m33 = 1 ;
308
-
309
- if (determinant (pmat) == 0 .f )
310
- {
311
- scale.x = 0 .f ;
312
- scale.y = 0 .f ;
313
- scale.z = 0 .f ;
314
-
315
- shear.x = 0 .f ;
316
- shear.y = 0 .f ;
317
- shear.z = 0 .f ;
318
-
319
- rotation.x = 0 .f ;
320
- rotation.y = 0 .f ;
321
- rotation.z = 0 .f ;
322
-
323
- perspective.x = perspective.y = perspective.z =
324
- perspective.w = 0 ;
325
-
326
- translation.x = locmat.m30 ;
327
- translation.y = locmat.m31 ;
328
- translation.z = locmat.m32 ;
329
- return 1 ;
330
- }
331
-
332
- /* First, isolate perspective. This is the messiest. */
333
- if (locmat.m03 != 0 || locmat.m13 != 0 ||
334
- locmat.m23 != 0 ) {
335
- /* prhs is the right hand side of the equation. */
336
- prhs.x = locmat.m03 ;
337
- prhs.y = locmat.m13 ;
338
- prhs.z = locmat.m23 ;
339
- prhs.w = locmat.m33 ;
340
-
341
- /* Solve the equation by inverting pmat and multiplying
342
- * prhs by the inverse. (This is the easiest way, not
343
- * necessarily the best.)
344
- * inverse function (and det4x4, above) from the Matrix
345
- * Inversion gem in the first volume.
346
- */
347
- invpmat = inverse (pmat);
348
- tinvpmat = invpmat.transpose ();
349
- psol = tinvpmat * prhs;
350
-
351
- /* Stuff the answer away. */
352
- perspective.x = psol.x ;
353
- perspective.y = psol.y ;
354
- perspective.z = psol.z ;
355
- perspective.w = psol.w ;
356
- /* Clear the perspective partition. */
357
- locmat.m03 = locmat.m13 =
358
- locmat.m23 = 0 ;
359
- locmat.m33 = 1 ;
360
- }
361
- else /* No perspective. */
362
- perspective.x = perspective.y = perspective.z =
363
- perspective.w = 0 ;
364
-
365
- /* Next take care of translation (easy). */
366
- translation.x = locmat.m30 ;
367
- translation.y = locmat.m31 ;
368
- translation.z = locmat.m32 ;
369
- locmat.m30 = 0 ;
370
- locmat.m31 = 0 ;
371
- locmat.m32 = 0 ;
372
-
373
-
374
- /* Now get scale and shear. */
375
- for (i = 0 ; i<3 ; i++) {
376
- row[i].x = locmat.m [i][0 ];
377
- row[i].y = locmat.m [i][1 ];
378
- row[i].z = locmat.m [i][2 ];
379
- }
380
-
381
- /* Compute X scale factor and normalize first row. */
382
- scale.x = sqrtf (row[0 ].sqnorm ());
383
- row[0 ].normalize ();
384
-
385
- /* Compute XY shear factor and make 2nd row orthogonal to 1st. */
386
- shear.x = dot (row[0 ], row[1 ]);
387
- row[1 ] = row[1 ] - shear.x * row[0 ];
388
-
389
- /* Now, compute Y scale and normalize 2nd row. */
390
- scale.y = sqrtf (row[1 ].sqnorm ());
391
- row[1 ].normalize ();
392
- shear.x /= scale.y ;
393
-
394
- /* Compute XZ and YZ shears, orthogonalize 3rd row. */
395
- shear.y = dot (row[0 ], row[2 ]);
396
- row[2 ] = row[2 ] - shear.y * row[0 ];
397
-
398
- shear.z = dot (row[1 ], row[2 ]);
399
- row[2 ] = row[2 ] - shear.z * row[1 ];
400
-
401
- /* Next, get Z scale and normalize 3rd row. */
402
- scale.z = sqrtf (row[2 ].sqnorm ());
403
- row[2 ].normalize ();
404
- shear.y /= scale.z ;
405
- shear.z /= scale.z ;
406
-
407
- /* At this point, the matrix (in rows[]) is orthonormal.
408
- * Check for a coordinate system flip. If the determinant
409
- * is -1, then negate the matrix and the scaling factors.
410
- */
411
- float3 r0 = float3 (row[0 ].x , row[0 ].y , row[0 ].z );
412
- float3 r1 = float3 (row[1 ].x , row[1 ].y , row[1 ].z );
413
- float3 r2 = float3 (row[2 ].x , row[2 ].y , row[2 ].z );
414
- if (dot (r0, cross (r1, r2)) < 0 )
415
- {
416
- for (i = 0 ; i < 3 ; i++) {
417
- row[i].x *= -1 ;
418
- row[i].y *= -1 ;
419
- row[i].z *= -1 ;
420
- }
421
- scale *= -1 ;
422
- }
423
-
424
-
425
- // source : http://www.gregslabaugh.net/publications/euler.pdf
426
- const float delta = 0 .000001f ;
427
- if ( fabsf (row[0 ].z ) < 1.0 - delta ) // if row[0].z is not +1 or -1
428
- {
429
- rotation.y = asinf (-row[0 ].z );
430
- const float cos_roty = cosf (rotation.y );
431
- rotation.x = atan2f (row[1 ].z / cos_roty, row[2 ].z / cos_roty);
432
- rotation.z = atan2f (row[0 ].y / cos_roty, row[0 ].x / cos_roty);
433
- }
434
- else
435
- {
436
- rotation.z = 0 ;
437
- if ( row[0 ].z < - (1.0 - delta) ) // if row[0].z =~ -1
438
- {
439
- rotation.y = asinf (-row[0 ].z ); // =~ + pi / 2.0;
440
- rotation.x = atan2f (row[1 ].x , row[2 ].x );
441
- }
442
- else
443
- {
444
- rotation.y = asinf (-row[0 ].z ); // =~ - pi / 2.0
445
- rotation.x = atan2f (-row[1 ].x , -row[2 ].x );
446
- }
447
- }
448
-
449
-
450
- /* All done! */
451
- return 1 ;
452
- }
453
275
}
0 commit comments