git.oblomov.eu Git - wine/blob - dlls/d3dx9_36/math.c

   1 /*
   2  * Mathematical operations specific to D3DX9.
   3  *
   4  * Copyright (C) 2008 David Adam
   5  * Copyright (C) 2008 Luis Busquets
   6  * Copyright (C) 2008 Jérôme Gardou
   7  * Copyright (C) 2008 Philip Nilsson
   8  * Copyright (C) 2008 Henri Verbeet
   9  *
  10  * This library is free software; you can redistribute it and/or
  11  * modify it under the terms of the GNU Lesser General Public
  12  * License as published by the Free Software Foundation; either
  13  * version 2.1 of the License, or (at your option) any later version.
  14  *
  15  * This library is distributed in the hope that it will be useful,
  16  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  17  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  18  * Lesser General Public License for more details.
  19  *
  20  * You should have received a copy of the GNU Lesser General Public
  21  * License along with this library; if not, write to the Free Software
  22  * Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA
  23  */
  24
  25 #define NONAMELESSUNION
  26
  27 #include "config.h"
  28 #include "wine/port.h"
  29
  30 #include "windef.h"
  31 #include "wingdi.h"
  32 #include "d3dx9_36_private.h"
  33
  34 #include "wine/debug.h"
  35
  36 WINE_DEFAULT_DEBUG_CHANNEL(d3dx);
  37
  38 static const ID3DXMatrixStackVtbl ID3DXMatrixStack_Vtbl;
  39
  40 struct ID3DXMatrixStackImpl
  41 {
  42   ID3DXMatrixStack ID3DXMatrixStack_iface;
  43   LONG ref;
  44
  45   unsigned int current;
  46   unsigned int stack_size;
  47   D3DXMATRIX *stack;
  48 };
  49
  50
  51 /*_________________D3DXColor____________________*/
  52
  53 D3DXCOLOR* WINAPI D3DXColorAdjustContrast(D3DXCOLOR *pout, const D3DXCOLOR *pc, FLOAT s)
  54 {
  55     TRACE("pout %p, pc %p, s %f\n", pout, pc, s);
  56
  57     pout->r = 0.5f + s * (pc->r - 0.5f);
  58     pout->g = 0.5f + s * (pc->g - 0.5f);
  59     pout->b = 0.5f + s * (pc->b - 0.5f);
  60     pout->a = pc->a;
  61     return pout;
  62 }
  63
  64 D3DXCOLOR* WINAPI D3DXColorAdjustSaturation(D3DXCOLOR *pout, const D3DXCOLOR *pc, FLOAT s)
  65 {
  66     FLOAT grey;
  67
  68     TRACE("pout %p, pc %p, s %f\n", pout, pc, s);
  69
  70     grey = pc->r * 0.2125f + pc->g * 0.7154f + pc->b * 0.0721f;
  71     pout->r = grey + s * (pc->r - grey);
  72     pout->g = grey + s * (pc->g - grey);
  73     pout->b = grey + s * (pc->b - grey);
  74     pout->a = pc->a;
  75     return pout;
  76 }
  77
  78 /*_________________Misc__________________________*/
  79
  80 FLOAT WINAPI D3DXFresnelTerm(FLOAT costheta, FLOAT refractionindex)
  81 {
  82     FLOAT a, d, g, result;
  83
  84     TRACE("costheta %f, refractionindex %f\n", costheta, refractionindex);
  85
  86     g = sqrtf(refractionindex * refractionindex + costheta * costheta - 1.0f);
  87     a = g + costheta;
  88     d = g - costheta;
  89     result = (costheta * a - 1.0f) * (costheta * a - 1.0f) / ((costheta * d + 1.0f) * (costheta * d + 1.0f)) + 1.0f;
  90     result *= 0.5f * d * d / (a * a);
  91
  92     return result;
  93 }
  94
  95 /*_________________D3DXMatrix____________________*/
  96
  97 D3DXMATRIX * WINAPI D3DXMatrixAffineTransformation(D3DXMATRIX *out, FLOAT scaling, const D3DXVECTOR3 *rotationcenter,
  98         const D3DXQUATERNION *rotation, const D3DXVECTOR3 *translation)
  99 {
 100     TRACE("out %p, scaling %f, rotationcenter %p, rotation %p, translation %p\n",
 101             out, scaling, rotationcenter, rotation, translation);
 102
 103     D3DXMatrixIdentity(out);
 104
 105     if (rotationcenter)
 106     {
 107         out->u.m[3][0] = -rotationcenter->x;
 108         out->u.m[3][1] = -rotationcenter->y;
 109         out->u.m[3][2] = -rotationcenter->z;
 110     }
 111
 112     if (rotation)
 113     {
 114         FLOAT temp00, temp01, temp02, temp10, temp11, temp12, temp20, temp21, temp22;
 115
 116         temp00 = 1.0f - 2.0f * (rotation->y * rotation->y + rotation->z * rotation->z);
 117         temp01 = 2.0f * (rotation->x * rotation->y + rotation->z * rotation->w);
 118         temp02 = 2.0f * (rotation->x * rotation->z - rotation->y * rotation->w);
 119         temp10 = 2.0f * (rotation->x * rotation->y - rotation->z * rotation->w);
 120         temp11 = 1.0f - 2.0f * (rotation->x * rotation->x + rotation->z * rotation->z);
 121         temp12 = 2.0f * (rotation->y * rotation->z + rotation->x * rotation->w);
 122         temp20 = 2.0f * (rotation->x * rotation->z + rotation->y * rotation->w);
 123         temp21 = 2.0f * (rotation->y * rotation->z - rotation->x * rotation->w);
 124         temp22 = 1.0f - 2.0f * (rotation->x * rotation->x + rotation->y * rotation->y);
 125
 126         out->u.m[0][0] = scaling * temp00;
 127         out->u.m[0][1] = scaling * temp01;
 128         out->u.m[0][2] = scaling * temp02;
 129         out->u.m[1][0] = scaling * temp10;
 130         out->u.m[1][1] = scaling * temp11;
 131         out->u.m[1][2] = scaling * temp12;
 132         out->u.m[2][0] = scaling * temp20;
 133         out->u.m[2][1] = scaling * temp21;
 134         out->u.m[2][2] = scaling * temp22;
 135
 136         if (rotationcenter)
 137         {
 138             FLOAT x, y, z;
 139
 140             x = out->u.m[3][0];
 141             y = out->u.m[3][1];
 142             z = out->u.m[3][2];
 143
 144             out->u.m[3][0] = x * temp00 + y * temp10 + z * temp20;
 145             out->u.m[3][1] = x * temp01 + y * temp11 + z * temp21;
 146             out->u.m[3][2] = x * temp02 + y * temp12 + z * temp22;
 147         }
 148     }
 149     else
 150     {
 151         out->u.m[0][0] = scaling;
 152         out->u.m[1][1] = scaling;
 153         out->u.m[2][2] = scaling;
 154     }
 155
 156     if (rotationcenter)
 157     {
 158         out->u.m[3][0] += rotationcenter->x;
 159         out->u.m[3][1] += rotationcenter->y;
 160         out->u.m[3][2] += rotationcenter->z;
 161     }
 162
 163     if (translation)
 164     {
 165         out->u.m[3][0] += translation->x;
 166         out->u.m[3][1] += translation->y;
 167         out->u.m[3][2] += translation->z;
 168     }
 169
 170     return out;
 171 }
 172
 173 D3DXMATRIX * WINAPI D3DXMatrixAffineTransformation2D(D3DXMATRIX *out, FLOAT scaling,
 174         const D3DXVECTOR2 *rotationcenter, FLOAT rotation, const D3DXVECTOR2 *translation)
 175 {
 176     FLOAT tmp1, tmp2, s;
 177
 178     TRACE("out %p, scaling %f, rotationcenter %p, rotation %f, translation %p\n",
 179             out, scaling, rotationcenter, rotation, translation);
 180
 181     s = sinf(rotation / 2.0f);
 182     tmp1 = 1.0f - 2.0f * s * s;
 183     tmp2 = 2.0 * s * cosf(rotation / 2.0f);
 184
 185     D3DXMatrixIdentity(out);
 186     out->u.m[0][0] = scaling * tmp1;
 187     out->u.m[0][1] = scaling * tmp2;
 188     out->u.m[1][0] = -scaling * tmp2;
 189     out->u.m[1][1] = scaling * tmp1;
 190
 191     if (rotationcenter)
 192     {
 193         FLOAT x, y;
 194
 195         x = rotationcenter->x;
 196         y = rotationcenter->y;
 197
 198         out->u.m[3][0] = y * tmp2 - x * tmp1 + x;
 199         out->u.m[3][1] = -x * tmp2 - y * tmp1 + y;
 200     }
 201
 202     if (translation)
 203     {
 204         out->u.m[3][0] += translation->x;
 205         out->u.m[3][1] += translation->y;
 206     }
 207
 208     return out;
 209 }
 210
 211 HRESULT WINAPI D3DXMatrixDecompose(D3DXVECTOR3 *poutscale, D3DXQUATERNION *poutrotation, D3DXVECTOR3 *pouttranslation, const D3DXMATRIX *pm)
 212 {
 213     D3DXMATRIX normalized;
 214     D3DXVECTOR3 vec;
 215
 216     TRACE("poutscale %p, poutrotation %p, pouttranslation %p, pm %p\n", poutscale, poutrotation, pouttranslation, pm);
 217
 218     /*Compute the scaling part.*/
 219     vec.x=pm->u.m[0][0];
 220     vec.y=pm->u.m[0][1];
 221     vec.z=pm->u.m[0][2];
 222     poutscale->x=D3DXVec3Length(&vec);
 223
 224     vec.x=pm->u.m[1][0];
 225     vec.y=pm->u.m[1][1];
 226     vec.z=pm->u.m[1][2];
 227     poutscale->y=D3DXVec3Length(&vec);
 228
 229     vec.x=pm->u.m[2][0];
 230     vec.y=pm->u.m[2][1];
 231     vec.z=pm->u.m[2][2];
 232     poutscale->z=D3DXVec3Length(&vec);
 233
 234     /*Compute the translation part.*/
 235     pouttranslation->x=pm->u.m[3][0];
 236     pouttranslation->y=pm->u.m[3][1];
 237     pouttranslation->z=pm->u.m[3][2];
 238
 239     /*Let's calculate the rotation now*/
 240     if ( (poutscale->x == 0.0f) || (poutscale->y == 0.0f) || (poutscale->z == 0.0f) ) return D3DERR_INVALIDCALL;
 241
 242     normalized.u.m[0][0]=pm->u.m[0][0]/poutscale->x;
 243     normalized.u.m[0][1]=pm->u.m[0][1]/poutscale->x;
 244     normalized.u.m[0][2]=pm->u.m[0][2]/poutscale->x;
 245     normalized.u.m[1][0]=pm->u.m[1][0]/poutscale->y;
 246     normalized.u.m[1][1]=pm->u.m[1][1]/poutscale->y;
 247     normalized.u.m[1][2]=pm->u.m[1][2]/poutscale->y;
 248     normalized.u.m[2][0]=pm->u.m[2][0]/poutscale->z;
 249     normalized.u.m[2][1]=pm->u.m[2][1]/poutscale->z;
 250     normalized.u.m[2][2]=pm->u.m[2][2]/poutscale->z;
 251
 252     D3DXQuaternionRotationMatrix(poutrotation,&normalized);
 253     return S_OK;
 254 }
 255
 256 FLOAT WINAPI D3DXMatrixDeterminant(const D3DXMATRIX *pm)
 257 {
 258     FLOAT t[3], v[4];
 259
 260     TRACE("pm %p\n", pm);
 261
 262     t[0] = pm->u.m[2][2] * pm->u.m[3][3] - pm->u.m[2][3] * pm->u.m[3][2];
 263     t[1] = pm->u.m[1][2] * pm->u.m[3][3] - pm->u.m[1][3] * pm->u.m[3][2];
 264     t[2] = pm->u.m[1][2] * pm->u.m[2][3] - pm->u.m[1][3] * pm->u.m[2][2];
 265     v[0] = pm->u.m[1][1] * t[0] - pm->u.m[2][1] * t[1] + pm->u.m[3][1] * t[2];
 266     v[1] = -pm->u.m[1][0] * t[0] + pm->u.m[2][0] * t[1] - pm->u.m[3][0] * t[2];
 267
 268     t[0] = pm->u.m[1][0] * pm->u.m[2][1] - pm->u.m[2][0] * pm->u.m[1][1];
 269     t[1] = pm->u.m[1][0] * pm->u.m[3][1] - pm->u.m[3][0] * pm->u.m[1][1];
 270     t[2] = pm->u.m[2][0] * pm->u.m[3][1] - pm->u.m[3][0] * pm->u.m[2][1];
 271     v[2] = pm->u.m[3][3] * t[0] - pm->u.m[2][3] * t[1] + pm->u.m[1][3] * t[2];
 272     v[3] = -pm->u.m[3][2] * t[0] + pm->u.m[2][2] * t[1] - pm->u.m[1][2] * t[2];
 273
 274     return pm->u.m[0][0] * v[0] + pm->u.m[0][1] * v[1] +
 275         pm->u.m[0][2] * v[2] + pm->u.m[0][3] * v[3];
 276 }
 277
 278 D3DXMATRIX* WINAPI D3DXMatrixInverse(D3DXMATRIX *pout, FLOAT *pdeterminant, const D3DXMATRIX *pm)
 279 {
 280     FLOAT det, t[3], v[16];
 281     UINT i, j;
 282
 283     TRACE("pout %p, pdeterminant %p, pm %p\n", pout, pdeterminant, pm);
 284
 285     t[0] = pm->u.m[2][2] * pm->u.m[3][3] - pm->u.m[2][3] * pm->u.m[3][2];
 286     t[1] = pm->u.m[1][2] * pm->u.m[3][3] - pm->u.m[1][3] * pm->u.m[3][2];
 287     t[2] = pm->u.m[1][2] * pm->u.m[2][3] - pm->u.m[1][3] * pm->u.m[2][2];
 288     v[0] = pm->u.m[1][1] * t[0] - pm->u.m[2][1] * t[1] + pm->u.m[3][1] * t[2];
 289     v[4] = -pm->u.m[1][0] * t[0] + pm->u.m[2][0] * t[1] - pm->u.m[3][0] * t[2];
 290
 291     t[0] = pm->u.m[1][0] * pm->u.m[2][1] - pm->u.m[2][0] * pm->u.m[1][1];
 292     t[1] = pm->u.m[1][0] * pm->u.m[3][1] - pm->u.m[3][0] * pm->u.m[1][1];
 293     t[2] = pm->u.m[2][0] * pm->u.m[3][1] - pm->u.m[3][0] * pm->u.m[2][1];
 294     v[8] = pm->u.m[3][3] * t[0] - pm->u.m[2][3] * t[1] + pm->u.m[1][3] * t[2];
 295     v[12] = -pm->u.m[3][2] * t[0] + pm->u.m[2][2] * t[1] - pm->u.m[1][2] * t[2];
 296
 297     det = pm->u.m[0][0] * v[0] + pm->u.m[0][1] * v[4] +
 298         pm->u.m[0][2] * v[8] + pm->u.m[0][3] * v[12];
 299     if (det == 0.0f)
 300         return NULL;
 301     if (pdeterminant)
 302         *pdeterminant = det;
 303
 304     t[0] = pm->u.m[2][2] * pm->u.m[3][3] - pm->u.m[2][3] * pm->u.m[3][2];
 305     t[1] = pm->u.m[0][2] * pm->u.m[3][3] - pm->u.m[0][3] * pm->u.m[3][2];
 306     t[2] = pm->u.m[0][2] * pm->u.m[2][3] - pm->u.m[0][3] * pm->u.m[2][2];
 307     v[1] = -pm->u.m[0][1] * t[0] + pm->u.m[2][1] * t[1] - pm->u.m[3][1] * t[2];
 308     v[5] = pm->u.m[0][0] * t[0] - pm->u.m[2][0] * t[1] + pm->u.m[3][0] * t[2];
 309
 310     t[0] = pm->u.m[0][0] * pm->u.m[2][1] - pm->u.m[2][0] * pm->u.m[0][1];
 311     t[1] = pm->u.m[3][0] * pm->u.m[0][1] - pm->u.m[0][0] * pm->u.m[3][1];
 312     t[2] = pm->u.m[2][0] * pm->u.m[3][1] - pm->u.m[3][0] * pm->u.m[2][1];
 313     v[9] = -pm->u.m[3][3] * t[0] - pm->u.m[2][3] * t[1]- pm->u.m[0][3] * t[2];
 314     v[13] = pm->u.m[3][2] * t[0] + pm->u.m[2][2] * t[1] + pm->u.m[0][2] * t[2];
 315
 316     t[0] = pm->u.m[1][2] * pm->u.m[3][3] - pm->u.m[1][3] * pm->u.m[3][2];
 317     t[1] = pm->u.m[0][2] * pm->u.m[3][3] - pm->u.m[0][3] * pm->u.m[3][2];
 318     t[2] = pm->u.m[0][2] * pm->u.m[1][3] - pm->u.m[0][3] * pm->u.m[1][2];
 319     v[2] = pm->u.m[0][1] * t[0] - pm->u.m[1][1] * t[1] + pm->u.m[3][1] * t[2];
 320     v[6] = -pm->u.m[0][0] * t[0] + pm->u.m[1][0] * t[1] - pm->u.m[3][0] * t[2];
 321
 322     t[0] = pm->u.m[0][0] * pm->u.m[1][1] - pm->u.m[1][0] * pm->u.m[0][1];
 323     t[1] = pm->u.m[3][0] * pm->u.m[0][1] - pm->u.m[0][0] * pm->u.m[3][1];
 324     t[2] = pm->u.m[1][0] * pm->u.m[3][1] - pm->u.m[3][0] * pm->u.m[1][1];
 325     v[10] = pm->u.m[3][3] * t[0] + pm->u.m[1][3] * t[1] + pm->u.m[0][3] * t[2];
 326     v[14] = -pm->u.m[3][2] * t[0] - pm->u.m[1][2] * t[1] - pm->u.m[0][2] * t[2];
 327
 328     t[0] = pm->u.m[1][2] * pm->u.m[2][3] - pm->u.m[1][3] * pm->u.m[2][2];
 329     t[1] = pm->u.m[0][2] * pm->u.m[2][3] - pm->u.m[0][3] * pm->u.m[2][2];
 330     t[2] = pm->u.m[0][2] * pm->u.m[1][3] - pm->u.m[0][3] * pm->u.m[1][2];
 331     v[3] = -pm->u.m[0][1] * t[0] + pm->u.m[1][1] * t[1] - pm->u.m[2][1] * t[2];
 332     v[7] = pm->u.m[0][0] * t[0] - pm->u.m[1][0] * t[1] + pm->u.m[2][0] * t[2];
 333
 334     v[11] = -pm->u.m[0][0] * (pm->u.m[1][1] * pm->u.m[2][3] - pm->u.m[1][3] * pm->u.m[2][1]) +
 335         pm->u.m[1][0] * (pm->u.m[0][1] * pm->u.m[2][3] - pm->u.m[0][3] * pm->u.m[2][1]) -
 336         pm->u.m[2][0] * (pm->u.m[0][1] * pm->u.m[1][3] - pm->u.m[0][3] * pm->u.m[1][1]);
 337
 338     v[15] = pm->u.m[0][0] * (pm->u.m[1][1] * pm->u.m[2][2] - pm->u.m[1][2] * pm->u.m[2][1]) -
 339         pm->u.m[1][0] * (pm->u.m[0][1] * pm->u.m[2][2] - pm->u.m[0][2] * pm->u.m[2][1]) +
 340         pm->u.m[2][0] * (pm->u.m[0][1] * pm->u.m[1][2] - pm->u.m[0][2] * pm->u.m[1][1]);
 341
 342     det = 1.0f / det;
 343
 344     for (i = 0; i < 4; i++)
 345         for (j = 0; j < 4; j++)
 346             pout->u.m[i][j] = v[4 * i + j] * det;
 347
 348     return pout;
 349 }
 350
 351 D3DXMATRIX* WINAPI D3DXMatrixLookAtLH(D3DXMATRIX *pout, const D3DXVECTOR3 *peye, const D3DXVECTOR3 *pat, const D3DXVECTOR3 *pup)
 352 {
 353     D3DXVECTOR3 right, rightn, up, upn, vec, vec2;
 354
 355     TRACE("pout %p, peye %p, pat %p, pup %p\n", pout, peye, pat, pup);
 356
 357     D3DXVec3Subtract(&vec2, pat, peye);
 358     D3DXVec3Normalize(&vec, &vec2);
 359     D3DXVec3Cross(&right, pup, &vec);
 360     D3DXVec3Cross(&up, &vec, &right);
 361     D3DXVec3Normalize(&rightn, &right);
 362     D3DXVec3Normalize(&upn, &up);
 363     pout->u.m[0][0] = rightn.x;
 364     pout->u.m[1][0] = rightn.y;
 365     pout->u.m[2][0] = rightn.z;
 366     pout->u.m[3][0] = -D3DXVec3Dot(&rightn,peye);
 367     pout->u.m[0][1] = upn.x;
 368     pout->u.m[1][1] = upn.y;
 369     pout->u.m[2][1] = upn.z;
 370     pout->u.m[3][1] = -D3DXVec3Dot(&upn, peye);
 371     pout->u.m[0][2] = vec.x;
 372     pout->u.m[1][2] = vec.y;
 373     pout->u.m[2][2] = vec.z;
 374     pout->u.m[3][2] = -D3DXVec3Dot(&vec, peye);
 375     pout->u.m[0][3] = 0.0f;
 376     pout->u.m[1][3] = 0.0f;
 377     pout->u.m[2][3] = 0.0f;
 378     pout->u.m[3][3] = 1.0f;
 379     return pout;
 380 }
 381
 382 D3DXMATRIX* WINAPI D3DXMatrixLookAtRH(D3DXMATRIX *pout, const D3DXVECTOR3 *peye, const D3DXVECTOR3 *pat, const D3DXVECTOR3 *pup)
 383 {
 384     D3DXVECTOR3 right, rightn, up, upn, vec, vec2;
 385
 386     TRACE("pout %p, peye %p, pat %p, pup %p\n", pout, peye, pat, pup);
 387
 388     D3DXVec3Subtract(&vec2, pat, peye);
 389     D3DXVec3Normalize(&vec, &vec2);
 390     D3DXVec3Cross(&right, pup, &vec);
 391     D3DXVec3Cross(&up, &vec, &right);
 392     D3DXVec3Normalize(&rightn, &right);
 393     D3DXVec3Normalize(&upn, &up);
 394     pout->u.m[0][0] = -rightn.x;
 395     pout->u.m[1][0] = -rightn.y;
 396     pout->u.m[2][0] = -rightn.z;
 397     pout->u.m[3][0] = D3DXVec3Dot(&rightn,peye);
 398     pout->u.m[0][1] = upn.x;
 399     pout->u.m[1][1] = upn.y;
 400     pout->u.m[2][1] = upn.z;
 401     pout->u.m[3][1] = -D3DXVec3Dot(&upn, peye);
 402     pout->u.m[0][2] = -vec.x;
 403     pout->u.m[1][2] = -vec.y;
 404     pout->u.m[2][2] = -vec.z;
 405     pout->u.m[3][2] = D3DXVec3Dot(&vec, peye);
 406     pout->u.m[0][3] = 0.0f;
 407     pout->u.m[1][3] = 0.0f;
 408     pout->u.m[2][3] = 0.0f;
 409     pout->u.m[3][3] = 1.0f;
 410     return pout;
 411 }
 412
 413 D3DXMATRIX* WINAPI D3DXMatrixMultiply(D3DXMATRIX *pout, const D3DXMATRIX *pm1, const D3DXMATRIX *pm2)
 414 {
 415     D3DXMATRIX out;
 416     int i,j;
 417
 418     TRACE("pout %p, pm1 %p, pm2 %p\n", pout, pm1, pm2);
 419
 420     for (i=0; i<4; i++)
 421     {
 422         for (j=0; j<4; j++)
 423         {
 424             out.u.m[i][j] = pm1->u.m[i][0] * pm2->u.m[0][j] + pm1->u.m[i][1] * pm2->u.m[1][j] + pm1->u.m[i][2] * pm2->u.m[2][j] + pm1->u.m[i][3] * pm2->u.m[3][j];
 425         }
 426     }
 427
 428     *pout = out;
 429     return pout;
 430 }
 431
 432 D3DXMATRIX* WINAPI D3DXMatrixMultiplyTranspose(D3DXMATRIX *pout, const D3DXMATRIX *pm1, const D3DXMATRIX *pm2)
 433 {
 434     D3DXMATRIX temp;
 435     int i, j;
 436
 437     TRACE("pout %p, pm1 %p, pm2 %p\n", pout, pm1, pm2);
 438
 439     for (i = 0; i < 4; i++)
 440         for (j = 0; j < 4; j++)
 441             temp.u.m[j][i] = pm1->u.m[i][0] * pm2->u.m[0][j] + pm1->u.m[i][1] * pm2->u.m[1][j] + pm1->u.m[i][2] * pm2->u.m[2][j] + pm1->u.m[i][3] * pm2->u.m[3][j];
 442
 443     *pout = temp;
 444     return pout;
 445 }
 446
 447 D3DXMATRIX* WINAPI D3DXMatrixOrthoLH(D3DXMATRIX *pout, FLOAT w, FLOAT h, FLOAT zn, FLOAT zf)
 448 {
 449     TRACE("pout %p, w %f, h %f, zn %f, zf %f\n", pout, w, h, zn, zf);
 450
 451     D3DXMatrixIdentity(pout);
 452     pout->u.m[0][0] = 2.0f / w;
 453     pout->u.m[1][1] = 2.0f / h;
 454     pout->u.m[2][2] = 1.0f / (zf - zn);
 455     pout->u.m[3][2] = zn / (zn - zf);
 456     return pout;
 457 }
 458
 459 D3DXMATRIX* WINAPI D3DXMatrixOrthoOffCenterLH(D3DXMATRIX *pout, FLOAT l, FLOAT r, FLOAT b, FLOAT t, FLOAT zn, FLOAT zf)
 460 {
 461     TRACE("pout %p, l %f, r %f, b %f, t %f, zn %f, zf %f\n", pout, l, r, b, t, zn, zf);
 462
 463     D3DXMatrixIdentity(pout);
 464     pout->u.m[0][0] = 2.0f / (r - l);
 465     pout->u.m[1][1] = 2.0f / (t - b);
 466     pout->u.m[2][2] = 1.0f / (zf -zn);
 467     pout->u.m[3][0] = -1.0f -2.0f *l / (r - l);
 468     pout->u.m[3][1] = 1.0f + 2.0f * t / (b - t);
 469     pout->u.m[3][2] = zn / (zn -zf);
 470     return pout;
 471 }
 472
 473 D3DXMATRIX* WINAPI D3DXMatrixOrthoOffCenterRH(D3DXMATRIX *pout, FLOAT l, FLOAT r, FLOAT b, FLOAT t, FLOAT zn, FLOAT zf)
 474 {
 475     TRACE("pout %p, l %f, r %f, b %f, t %f, zn %f, zf %f\n", pout, l, r, b, t, zn, zf);
 476
 477     D3DXMatrixIdentity(pout);
 478     pout->u.m[0][0] = 2.0f / (r - l);
 479     pout->u.m[1][1] = 2.0f / (t - b);
 480     pout->u.m[2][2] = 1.0f / (zn -zf);
 481     pout->u.m[3][0] = -1.0f -2.0f *l / (r - l);
 482     pout->u.m[3][1] = 1.0f + 2.0f * t / (b - t);
 483     pout->u.m[3][2] = zn / (zn -zf);
 484     return pout;
 485 }
 486
 487 D3DXMATRIX* WINAPI D3DXMatrixOrthoRH(D3DXMATRIX *pout, FLOAT w, FLOAT h, FLOAT zn, FLOAT zf)
 488 {
 489     TRACE("pout %p, w %f, h %f, zn %f, zf %f\n", pout, w, h, zn, zf);
 490
 491     D3DXMatrixIdentity(pout);
 492     pout->u.m[0][0] = 2.0f / w;
 493     pout->u.m[1][1] = 2.0f / h;
 494     pout->u.m[2][2] = 1.0f / (zn - zf);
 495     pout->u.m[3][2] = zn / (zn - zf);
 496     return pout;
 497 }
 498
 499 D3DXMATRIX* WINAPI D3DXMatrixPerspectiveFovLH(D3DXMATRIX *pout, FLOAT fovy, FLOAT aspect, FLOAT zn, FLOAT zf)
 500 {
 501     TRACE("pout %p, fovy %f, aspect %f, zn %f, zf %f\n", pout, fovy, aspect, zn, zf);
 502
 503     D3DXMatrixIdentity(pout);
 504     pout->u.m[0][0] = 1.0f / (aspect * tan(fovy/2.0f));
 505     pout->u.m[1][1] = 1.0f / tan(fovy/2.0f);
 506     pout->u.m[2][2] = zf / (zf - zn);
 507     pout->u.m[2][3] = 1.0f;
 508     pout->u.m[3][2] = (zf * zn) / (zn - zf);
 509     pout->u.m[3][3] = 0.0f;
 510     return pout;
 511 }
 512
 513 D3DXMATRIX* WINAPI D3DXMatrixPerspectiveFovRH(D3DXMATRIX *pout, FLOAT fovy, FLOAT aspect, FLOAT zn, FLOAT zf)
 514 {
 515     TRACE("pout %p, fovy %f, aspect %f, zn %f, zf %f\n", pout, fovy, aspect, zn, zf);
 516
 517     D3DXMatrixIdentity(pout);
 518     pout->u.m[0][0] = 1.0f / (aspect * tan(fovy/2.0f));
 519     pout->u.m[1][1] = 1.0f / tan(fovy/2.0f);
 520     pout->u.m[2][2] = zf / (zn - zf);
 521     pout->u.m[2][3] = -1.0f;
 522     pout->u.m[3][2] = (zf * zn) / (zn - zf);
 523     pout->u.m[3][3] = 0.0f;
 524     return pout;
 525 }
 526
 527 D3DXMATRIX* WINAPI D3DXMatrixPerspectiveLH(D3DXMATRIX *pout, FLOAT w, FLOAT h, FLOAT zn, FLOAT zf)
 528 {
 529     TRACE("pout %p, w %f, h %f, zn %f, zf %f\n", pout, w, h, zn, zf);
 530
 531     D3DXMatrixIdentity(pout);
 532     pout->u.m[0][0] = 2.0f * zn / w;
 533     pout->u.m[1][1] = 2.0f * zn / h;
 534     pout->u.m[2][2] = zf / (zf - zn);
 535     pout->u.m[3][2] = (zn * zf) / (zn - zf);
 536     pout->u.m[2][3] = 1.0f;
 537     pout->u.m[3][3] = 0.0f;
 538     return pout;
 539 }
 540
 541 D3DXMATRIX* WINAPI D3DXMatrixPerspectiveOffCenterLH(D3DXMATRIX *pout, FLOAT l, FLOAT r, FLOAT b, FLOAT t, FLOAT zn, FLOAT zf)
 542 {
 543     TRACE("pout %p, l %f, r %f, b %f, t %f, zn %f, zf %f\n", pout, l, r, b, t, zn, zf);
 544
 545     D3DXMatrixIdentity(pout);
 546     pout->u.m[0][0] = 2.0f * zn / (r - l);
 547     pout->u.m[1][1] = -2.0f * zn / (b - t);
 548     pout->u.m[2][0] = -1.0f - 2.0f * l / (r - l);
 549     pout->u.m[2][1] = 1.0f + 2.0f * t / (b - t);
 550     pout->u.m[2][2] = - zf / (zn - zf);
 551     pout->u.m[3][2] = (zn * zf) / (zn -zf);
 552     pout->u.m[2][3] = 1.0f;
 553     pout->u.m[3][3] = 0.0f;
 554     return pout;
 555 }
 556
 557 D3DXMATRIX* WINAPI D3DXMatrixPerspectiveOffCenterRH(D3DXMATRIX *pout, FLOAT l, FLOAT r, FLOAT b, FLOAT t, FLOAT zn, FLOAT zf)
 558 {
 559     TRACE("pout %p, l %f, r %f, b %f, t %f, zn %f, zf %f\n", pout, l, r, b, t, zn, zf);
 560
 561     D3DXMatrixIdentity(pout);
 562     pout->u.m[0][0] = 2.0f * zn / (r - l);
 563     pout->u.m[1][1] = -2.0f * zn / (b - t);
 564     pout->u.m[2][0] = 1.0f + 2.0f * l / (r - l);
 565     pout->u.m[2][1] = -1.0f -2.0f * t / (b - t);
 566     pout->u.m[2][2] = zf / (zn - zf);
 567     pout->u.m[3][2] = (zn * zf) / (zn -zf);
 568     pout->u.m[2][3] = -1.0f;
 569     pout->u.m[3][3] = 0.0f;
 570     return pout;
 571 }
 572
 573 D3DXMATRIX* WINAPI D3DXMatrixPerspectiveRH(D3DXMATRIX *pout, FLOAT w, FLOAT h, FLOAT zn, FLOAT zf)
 574 {
 575     TRACE("pout %p, w %f, h %f, zn %f, zf %f\n", pout, w, h, zn, zf);
 576
 577     D3DXMatrixIdentity(pout);
 578     pout->u.m[0][0] = 2.0f * zn / w;
 579     pout->u.m[1][1] = 2.0f * zn / h;
 580     pout->u.m[2][2] = zf / (zn - zf);
 581     pout->u.m[3][2] = (zn * zf) / (zn - zf);
 582     pout->u.m[2][3] = -1.0f;
 583     pout->u.m[3][3] = 0.0f;
 584     return pout;
 585 }
 586
 587 D3DXMATRIX* WINAPI D3DXMatrixReflect(D3DXMATRIX *pout, const D3DXPLANE *pplane)
 588 {
 589     D3DXPLANE Nplane;
 590
 591     TRACE("pout %p, pplane %p\n", pout, pplane);
 592
 593     D3DXPlaneNormalize(&Nplane, pplane);
 594     D3DXMatrixIdentity(pout);
 595     pout->u.m[0][0] = 1.0f - 2.0f * Nplane.a * Nplane.a;
 596     pout->u.m[0][1] = -2.0f * Nplane.a * Nplane.b;
 597     pout->u.m[0][2] = -2.0f * Nplane.a * Nplane.c;
 598     pout->u.m[1][0] = -2.0f * Nplane.a * Nplane.b;
 599     pout->u.m[1][1] = 1.0f - 2.0f * Nplane.b * Nplane.b;
 600     pout->u.m[1][2] = -2.0f * Nplane.b * Nplane.c;
 601     pout->u.m[2][0] = -2.0f * Nplane.c * Nplane.a;
 602     pout->u.m[2][1] = -2.0f * Nplane.c * Nplane.b;
 603     pout->u.m[2][2] = 1.0f - 2.0f * Nplane.c * Nplane.c;
 604     pout->u.m[3][0] = -2.0f * Nplane.d * Nplane.a;
 605     pout->u.m[3][1] = -2.0f * Nplane.d * Nplane.b;
 606     pout->u.m[3][2] = -2.0f * Nplane.d * Nplane.c;
 607     return pout;
 608 }
 609
 610 D3DXMATRIX * WINAPI D3DXMatrixRotationAxis(D3DXMATRIX *out, const D3DXVECTOR3 *v, FLOAT angle)
 611 {
 612     D3DXVECTOR3 nv;
 613     FLOAT sangle, cangle, cdiff;
 614
 615     TRACE("out %p, v %p, angle %f\n", out, v, angle);
 616
 617     D3DXVec3Normalize(&nv, v);
 618     sangle = sinf(angle);
 619     cangle = cosf(angle);
 620     cdiff = 1.0f - cangle;
 621
 622     out->u.m[0][0] = cdiff * nv.x * nv.x + cangle;
 623     out->u.m[1][0] = cdiff * nv.x * nv.y - sangle * nv.z;
 624     out->u.m[2][0] = cdiff * nv.x * nv.z + sangle * nv.y;
 625     out->u.m[3][0] = 0.0f;
 626     out->u.m[0][1] = cdiff * nv.y * nv.x + sangle * nv.z;
 627     out->u.m[1][1] = cdiff * nv.y * nv.y + cangle;
 628     out->u.m[2][1] = cdiff * nv.y * nv.z - sangle * nv.x;
 629     out->u.m[3][1] = 0.0f;
 630     out->u.m[0][2] = cdiff * nv.z * nv.x - sangle * nv.y;
 631     out->u.m[1][2] = cdiff * nv.z * nv.y + sangle * nv.x;
 632     out->u.m[2][2] = cdiff * nv.z * nv.z + cangle;
 633     out->u.m[3][2] = 0.0f;
 634     out->u.m[0][3] = 0.0f;
 635     out->u.m[1][3] = 0.0f;
 636     out->u.m[2][3] = 0.0f;
 637     out->u.m[3][3] = 1.0f;
 638
 639     return out;
 640 }
 641
 642 D3DXMATRIX* WINAPI D3DXMatrixRotationQuaternion(D3DXMATRIX *pout, const D3DXQUATERNION *pq)
 643 {
 644     TRACE("pout %p, pq %p\n", pout, pq);
 645
 646     D3DXMatrixIdentity(pout);
 647     pout->u.m[0][0] = 1.0f - 2.0f * (pq->y * pq->y + pq->z * pq->z);
 648     pout->u.m[0][1] = 2.0f * (pq->x *pq->y + pq->z * pq->w);
 649     pout->u.m[0][2] = 2.0f * (pq->x * pq->z - pq->y * pq->w);
 650     pout->u.m[1][0] = 2.0f * (pq->x * pq->y - pq->z * pq->w);
 651     pout->u.m[1][1] = 1.0f - 2.0f * (pq->x * pq->x + pq->z * pq->z);
 652     pout->u.m[1][2] = 2.0f * (pq->y *pq->z + pq->x *pq->w);
 653     pout->u.m[2][0] = 2.0f * (pq->x * pq->z + pq->y * pq->w);
 654     pout->u.m[2][1] = 2.0f * (pq->y *pq->z - pq->x *pq->w);
 655     pout->u.m[2][2] = 1.0f - 2.0f * (pq->x * pq->x + pq->y * pq->y);
 656     return pout;
 657 }
 658
 659 D3DXMATRIX* WINAPI D3DXMatrixRotationX(D3DXMATRIX *pout, FLOAT angle)
 660 {
 661     TRACE("pout %p, angle %f\n", pout, angle);
 662
 663     D3DXMatrixIdentity(pout);
 664     pout->u.m[1][1] = cos(angle);
 665     pout->u.m[2][2] = cos(angle);
 666     pout->u.m[1][2] = sin(angle);
 667     pout->u.m[2][1] = -sin(angle);
 668     return pout;
 669 }
 670
 671 D3DXMATRIX* WINAPI D3DXMatrixRotationY(D3DXMATRIX *pout, FLOAT angle)
 672 {
 673     TRACE("pout %p, angle %f\n", pout, angle);
 674
 675     D3DXMatrixIdentity(pout);
 676     pout->u.m[0][0] = cos(angle);
 677     pout->u.m[2][2] = cos(angle);
 678     pout->u.m[0][2] = -sin(angle);
 679     pout->u.m[2][0] = sin(angle);
 680     return pout;
 681 }
 682
 683 D3DXMATRIX * WINAPI D3DXMatrixRotationYawPitchRoll(D3DXMATRIX *out, FLOAT yaw, FLOAT pitch, FLOAT roll)
 684 {
 685     FLOAT sroll, croll, spitch, cpitch, syaw, cyaw;
 686
 687     TRACE("out %p, yaw %f, pitch %f, roll %f\n", out, yaw, pitch, roll);
 688
 689     sroll = sinf(roll);
 690     croll = cosf(roll);
 691     spitch = sinf(pitch);
 692     cpitch = cosf(pitch);
 693     syaw = sinf(yaw);
 694     cyaw = cosf(yaw);
 695
 696     out->u.m[0][0] = sroll * spitch * syaw + croll * cyaw;
 697     out->u.m[0][1] = sroll * cpitch;
 698     out->u.m[0][2] = sroll * spitch * cyaw - croll * syaw;
 699     out->u.m[0][3] = 0.0f;
 700     out->u.m[1][0] = croll * spitch * syaw - sroll * cyaw;
 701     out->u.m[1][1] = croll * cpitch;
 702     out->u.m[1][2] = croll * spitch * cyaw + sroll * syaw;
 703     out->u.m[1][3] = 0.0f;
 704     out->u.m[2][0] = cpitch * syaw;
 705     out->u.m[2][1] = -spitch;
 706     out->u.m[2][2] = cpitch * cyaw;
 707     out->u.m[2][3] = 0.0f;
 708     out->u.m[3][0] = 0.0f;
 709     out->u.m[3][1] = 0.0f;
 710     out->u.m[3][2] = 0.0f;
 711     out->u.m[3][3] = 1.0f;
 712
 713     return out;
 714 }
 715
 716 D3DXMATRIX* WINAPI D3DXMatrixRotationZ(D3DXMATRIX *pout, FLOAT angle)
 717 {
 718     TRACE("pout %p, angle %f\n", pout, angle);
 719
 720     D3DXMatrixIdentity(pout);
 721     pout->u.m[0][0] = cos(angle);
 722     pout->u.m[1][1] = cos(angle);
 723     pout->u.m[0][1] = sin(angle);
 724     pout->u.m[1][0] = -sin(angle);
 725     return pout;
 726 }
 727
 728 D3DXMATRIX* WINAPI D3DXMatrixScaling(D3DXMATRIX *pout, FLOAT sx, FLOAT sy, FLOAT sz)
 729 {
 730     TRACE("pout %p, sx %f, sy %f, sz %f\n", pout, sx, sy, sz);
 731
 732     D3DXMatrixIdentity(pout);
 733     pout->u.m[0][0] = sx;
 734     pout->u.m[1][1] = sy;
 735     pout->u.m[2][2] = sz;
 736     return pout;
 737 }
 738
 739 D3DXMATRIX* WINAPI D3DXMatrixShadow(D3DXMATRIX *pout, const D3DXVECTOR4 *plight, const D3DXPLANE *pplane)
 740 {
 741     D3DXPLANE Nplane;
 742     FLOAT dot;
 743
 744     TRACE("pout %p, plight %p, pplane %p\n", pout, plight, pplane);
 745
 746     D3DXPlaneNormalize(&Nplane, pplane);
 747     dot = D3DXPlaneDot(&Nplane, plight);
 748     pout->u.m[0][0] = dot - Nplane.a * plight->x;
 749     pout->u.m[0][1] = -Nplane.a * plight->y;
 750     pout->u.m[0][2] = -Nplane.a * plight->z;
 751     pout->u.m[0][3] = -Nplane.a * plight->w;
 752     pout->u.m[1][0] = -Nplane.b * plight->x;
 753     pout->u.m[1][1] = dot - Nplane.b * plight->y;
 754     pout->u.m[1][2] = -Nplane.b * plight->z;
 755     pout->u.m[1][3] = -Nplane.b * plight->w;
 756     pout->u.m[2][0] = -Nplane.c * plight->x;
 757     pout->u.m[2][1] = -Nplane.c * plight->y;
 758     pout->u.m[2][2] = dot - Nplane.c * plight->z;
 759     pout->u.m[2][3] = -Nplane.c * plight->w;
 760     pout->u.m[3][0] = -Nplane.d * plight->x;
 761     pout->u.m[3][1] = -Nplane.d * plight->y;
 762     pout->u.m[3][2] = -Nplane.d * plight->z;
 763     pout->u.m[3][3] = dot - Nplane.d * plight->w;
 764     return pout;
 765 }
 766
 767 D3DXMATRIX* WINAPI D3DXMatrixTransformation(D3DXMATRIX *pout, const D3DXVECTOR3 *pscalingcenter, const D3DXQUATERNION *pscalingrotation, const D3DXVECTOR3 *pscaling, const D3DXVECTOR3 *protationcenter, const D3DXQUATERNION *protation, const D3DXVECTOR3 *ptranslation)
 768 {
 769     D3DXMATRIX m1, m2, m3, m4, m5, m6, m7;
 770     D3DXQUATERNION prc;
 771     D3DXVECTOR3 psc, pt;
 772
 773     TRACE("pout %p, pscalingcenter %p, pscalingrotation %p, pscaling %p, protationcentr %p, protation %p, ptranslation %p\n",
 774         pout, pscalingcenter, pscalingrotation, pscaling, protationcenter, protation, ptranslation);
 775
 776     if ( !pscalingcenter )
 777     {
 778         psc.x = 0.0f;
 779         psc.y = 0.0f;
 780         psc.z = 0.0f;
 781     }
 782     else
 783     {
 784         psc.x = pscalingcenter->x;
 785         psc.y = pscalingcenter->y;
 786         psc.z = pscalingcenter->z;
 787     }
 788
 789     if ( !protationcenter )
 790     {
 791         prc.x = 0.0f;
 792         prc.y = 0.0f;
 793         prc.z = 0.0f;
 794     }
 795     else
 796     {
 797         prc.x = protationcenter->x;
 798         prc.y = protationcenter->y;
 799         prc.z = protationcenter->z;
 800     }
 801
 802     if ( !ptranslation )
 803     {
 804         pt.x = 0.0f;
 805         pt.y = 0.0f;
 806         pt.z = 0.0f;
 807     }
 808     else
 809     {
 810         pt.x = ptranslation->x;
 811         pt.y = ptranslation->y;
 812         pt.z = ptranslation->z;
 813     }
 814
 815     D3DXMatrixTranslation(&m1, -psc.x, -psc.y, -psc.z);
 816
 817     if ( !pscalingrotation )
 818     {
 819         D3DXMatrixIdentity(&m2);
 820         D3DXMatrixIdentity(&m4);
 821     }
 822     else
 823     {
 824         D3DXMatrixRotationQuaternion(&m4, pscalingrotation);
 825         D3DXMatrixInverse(&m2, NULL, &m4);
 826     }
 827
 828     if ( !pscaling ) D3DXMatrixIdentity(&m3);
 829     else D3DXMatrixScaling(&m3, pscaling->x, pscaling->y, pscaling->z);
 830
 831     if ( !protation ) D3DXMatrixIdentity(&m6);
 832     else D3DXMatrixRotationQuaternion(&m6, protation);
 833
 834     D3DXMatrixTranslation(&m5, psc.x - prc.x,  psc.y - prc.y,  psc.z - prc.z);
 835     D3DXMatrixTranslation(&m7, prc.x + pt.x, prc.y + pt.y, prc.z + pt.z);
 836     D3DXMatrixMultiply(&m1, &m1, &m2);
 837     D3DXMatrixMultiply(&m1, &m1, &m3);
 838     D3DXMatrixMultiply(&m1, &m1, &m4);
 839     D3DXMatrixMultiply(&m1, &m1, &m5);
 840     D3DXMatrixMultiply(&m1, &m1, &m6);
 841     D3DXMatrixMultiply(pout, &m1, &m7);
 842     return pout;
 843 }
 844
 845 D3DXMATRIX* WINAPI D3DXMatrixTransformation2D(D3DXMATRIX *pout, const D3DXVECTOR2 *pscalingcenter, FLOAT scalingrotation, const D3DXVECTOR2 *pscaling, const D3DXVECTOR2 *protationcenter, FLOAT rotation, const D3DXVECTOR2 *ptranslation)
 846 {
 847     D3DXQUATERNION rot, sca_rot;
 848     D3DXVECTOR3 rot_center, sca, sca_center, trans;
 849
 850     TRACE("pout %p, pscalingcenter %p, scalingrotation %f, pscaling %p, protztioncenter %p, rotation %f, ptranslation %p\n",
 851         pout, pscalingcenter, scalingrotation, pscaling, protationcenter, rotation, ptranslation);
 852
 853     if ( pscalingcenter )
 854     {
 855         sca_center.x=pscalingcenter->x;
 856         sca_center.y=pscalingcenter->y;
 857         sca_center.z=0.0f;
 858     }
 859     else
 860     {
 861         sca_center.x=0.0f;
 862         sca_center.y=0.0f;
 863         sca_center.z=0.0f;
 864     }
 865
 866     if ( pscaling )
 867     {
 868         sca.x=pscaling->x;
 869         sca.y=pscaling->y;
 870         sca.z=1.0f;
 871     }
 872     else
 873     {
 874         sca.x=1.0f;
 875         sca.y=1.0f;
 876         sca.z=1.0f;
 877     }
 878
 879     if ( protationcenter )
 880     {
 881         rot_center.x=protationcenter->x;
 882         rot_center.y=protationcenter->y;
 883         rot_center.z=0.0f;
 884     }
 885     else
 886     {
 887         rot_center.x=0.0f;
 888         rot_center.y=0.0f;
 889         rot_center.z=0.0f;
 890     }
 891
 892     if ( ptranslation )
 893     {
 894         trans.x=ptranslation->x;
 895         trans.y=ptranslation->y;
 896         trans.z=0.0f;
 897     }
 898     else
 899     {
 900         trans.x=0.0f;
 901         trans.y=0.0f;
 902         trans.z=0.0f;
 903     }
 904
 905     rot.w=cos(rotation/2.0f);
 906     rot.x=0.0f;
 907     rot.y=0.0f;
 908     rot.z=sin(rotation/2.0f);
 909
 910     sca_rot.w=cos(scalingrotation/2.0f);
 911     sca_rot.x=0.0f;
 912     sca_rot.y=0.0f;
 913     sca_rot.z=sin(scalingrotation/2.0f);
 914
 915     D3DXMatrixTransformation(pout, &sca_center, &sca_rot, &sca, &rot_center, &rot, &trans);
 916
 917     return pout;
 918 }
 919
 920 D3DXMATRIX* WINAPI D3DXMatrixTranslation(D3DXMATRIX *pout, FLOAT x, FLOAT y, FLOAT z)
 921 {
 922     TRACE("pout %p, x %f, y %f, z %f\n", pout, x, y, z);
 923
 924     D3DXMatrixIdentity(pout);
 925     pout->u.m[3][0] = x;
 926     pout->u.m[3][1] = y;
 927     pout->u.m[3][2] = z;
 928     return pout;
 929 }
 930
 931 D3DXMATRIX* WINAPI D3DXMatrixTranspose(D3DXMATRIX *pout, const D3DXMATRIX *pm)
 932 {
 933     const D3DXMATRIX m = *pm;
 934     int i,j;
 935
 936     TRACE("pout %p, pm %p\n", pout, pm);
 937
 938     for (i=0; i<4; i++)
 939         for (j=0; j<4; j++) pout->u.m[i][j] = m.u.m[j][i];
 940
 941     return pout;
 942 }
 943
 944 /*_________________D3DXMatrixStack____________________*/
 945
 946 static const unsigned int INITIAL_STACK_SIZE = 32;
 947
 948 HRESULT WINAPI D3DXCreateMatrixStack(DWORD flags, ID3DXMatrixStack **ppstack)
 949 {
 950     struct ID3DXMatrixStackImpl *object;
 951
 952     TRACE("flags %#x, ppstack %p\n", flags, ppstack);
 953
 954     object = HeapAlloc(GetProcessHeap(), HEAP_ZERO_MEMORY, sizeof(*object));
 955     if (object == NULL)
 956     {
 957         *ppstack = NULL;
 958         return E_OUTOFMEMORY;
 959     }
 960     object->ID3DXMatrixStack_iface.lpVtbl = &ID3DXMatrixStack_Vtbl;
 961     object->ref = 1;
 962
 963     object->stack = HeapAlloc(GetProcessHeap(), 0, INITIAL_STACK_SIZE * sizeof(*object->stack));
 964     if (!object->stack)
 965     {
 966         HeapFree(GetProcessHeap(), 0, object);
 967         *ppstack = NULL;
 968         return E_OUTOFMEMORY;
 969     }
 970
 971     object->current = 0;
 972     object->stack_size = INITIAL_STACK_SIZE;
 973     D3DXMatrixIdentity(&object->stack[0]);
 974
 975     TRACE("Created matrix stack %p\n", object);
 976
 977     *ppstack = &object->ID3DXMatrixStack_iface;
 978     return D3D_OK;
 979 }
 980
 981 static inline struct ID3DXMatrixStackImpl *impl_from_ID3DXMatrixStack(ID3DXMatrixStack *iface)
 982 {
 983   return CONTAINING_RECORD(iface, struct ID3DXMatrixStackImpl, ID3DXMatrixStack_iface);
 984 }
 985
 986 static HRESULT WINAPI ID3DXMatrixStackImpl_QueryInterface(ID3DXMatrixStack *iface, REFIID riid, void **out)
 987 {
 988     TRACE("iface %p, riid %s, out %p.\n", iface, debugstr_guid(riid), out);
 989
 990     if (IsEqualGUID(riid, &IID_ID3DXMatrixStack)
 991             || IsEqualGUID(riid, &IID_IUnknown))
 992     {
 993         ID3DXMatrixStack_AddRef(iface);
 994         *out = iface;
 995         return S_OK;
 996     }
 997
 998     WARN("%s not implemented, returning E_NOINTERFACE.\n", debugstr_guid(riid));
 999
1000     *out = NULL;
1001     return E_NOINTERFACE;
1002 }
1003
1004 static ULONG WINAPI ID3DXMatrixStackImpl_AddRef(ID3DXMatrixStack *iface)
1005 {
1006     struct ID3DXMatrixStackImpl *This = impl_from_ID3DXMatrixStack(iface);
1007     ULONG ref = InterlockedIncrement(&This->ref);
1008     TRACE("(%p) : AddRef from %d\n", This, ref - 1);
1009     return ref;
1010 }
1011
1012 static ULONG WINAPI ID3DXMatrixStackImpl_Release(ID3DXMatrixStack *iface)
1013 {
1014     struct ID3DXMatrixStackImpl *This = impl_from_ID3DXMatrixStack(iface);
1015     ULONG ref = InterlockedDecrement(&This->ref);
1016     if (!ref)
1017     {
1018         HeapFree(GetProcessHeap(), 0, This->stack);
1019         HeapFree(GetProcessHeap(), 0, This);
1020     }
1021     TRACE("(%p) : ReleaseRef to %d\n", This, ref);
1022     return ref;
1023 }
1024
1025 static D3DXMATRIX* WINAPI ID3DXMatrixStackImpl_GetTop(ID3DXMatrixStack *iface)
1026 {
1027     struct ID3DXMatrixStackImpl *This = impl_from_ID3DXMatrixStack(iface);
1028
1029     TRACE("iface %p\n", iface);
1030
1031     return &This->stack[This->current];
1032 }
1033
1034 static HRESULT WINAPI ID3DXMatrixStackImpl_LoadIdentity(ID3DXMatrixStack *iface)
1035 {
1036     struct ID3DXMatrixStackImpl *This = impl_from_ID3DXMatrixStack(iface);
1037
1038     TRACE("iface %p\n", iface);
1039
1040     D3DXMatrixIdentity(&This->stack[This->current]);
1041
1042     return D3D_OK;
1043 }
1044
1045 static HRESULT WINAPI ID3DXMatrixStackImpl_LoadMatrix(ID3DXMatrixStack *iface, const D3DXMATRIX *pm)
1046 {
1047     struct ID3DXMatrixStackImpl *This = impl_from_ID3DXMatrixStack(iface);
1048
1049     TRACE("iface %p, pm %p\n", iface, pm);
1050
1051     This->stack[This->current] = *pm;
1052
1053     return D3D_OK;
1054 }
1055
1056 static HRESULT WINAPI ID3DXMatrixStackImpl_MultMatrix(ID3DXMatrixStack *iface, const D3DXMATRIX *pm)
1057 {
1058     struct ID3DXMatrixStackImpl *This = impl_from_ID3DXMatrixStack(iface);
1059
1060     TRACE("iface %p, pm %p\n", iface, pm);
1061
1062     D3DXMatrixMultiply(&This->stack[This->current], &This->stack[This->current], pm);
1063
1064     return D3D_OK;
1065 }
1066
1067 static HRESULT WINAPI ID3DXMatrixStackImpl_MultMatrixLocal(ID3DXMatrixStack *iface, const D3DXMATRIX *pm)
1068 {
1069     struct ID3DXMatrixStackImpl *This = impl_from_ID3DXMatrixStack(iface);
1070
1071     TRACE("iface %p, pm %p\n", iface, pm);
1072
1073     D3DXMatrixMultiply(&This->stack[This->current], pm, &This->stack[This->current]);
1074
1075     return D3D_OK;
1076 }
1077
1078 static HRESULT WINAPI ID3DXMatrixStackImpl_Pop(ID3DXMatrixStack *iface)
1079 {
1080     struct ID3DXMatrixStackImpl *This = impl_from_ID3DXMatrixStack(iface);
1081
1082     TRACE("iface %p\n", iface);
1083
1084     /* Popping the last element on the stack returns D3D_OK, but does nothing. */
1085     if (!This->current) return D3D_OK;
1086
1087     if (This->current <= This->stack_size / 4 && This->stack_size >= INITIAL_STACK_SIZE * 2)
1088     {
1089         unsigned int new_size;
1090         D3DXMATRIX *new_stack;
1091
1092         new_size = This->stack_size / 2;
1093         new_stack = HeapReAlloc(GetProcessHeap(), 0, This->stack, new_size * sizeof(*new_stack));
1094         if (new_stack)
1095         {
1096             This->stack_size = new_size;
1097             This->stack = new_stack;
1098         }
1099     }
1100
1101     --This->current;
1102
1103     return D3D_OK;
1104 }
1105
1106 static HRESULT WINAPI ID3DXMatrixStackImpl_Push(ID3DXMatrixStack *iface)
1107 {
1108     struct ID3DXMatrixStackImpl *This = impl_from_ID3DXMatrixStack(iface);
1109
1110     TRACE("iface %p\n", iface);
1111
1112     if (This->current == This->stack_size - 1)
1113     {
1114         unsigned int new_size;
1115         D3DXMATRIX *new_stack;
1116
1117         if (This->stack_size > UINT_MAX / 2) return E_OUTOFMEMORY;
1118
1119         new_size = This->stack_size * 2;
1120         new_stack = HeapReAlloc(GetProcessHeap(), 0, This->stack, new_size * sizeof(*new_stack));
1121         if (!new_stack) return E_OUTOFMEMORY;
1122
1123         This->stack_size = new_size;
1124         This->stack = new_stack;
1125     }
1126
1127     ++This->current;
1128     This->stack[This->current] = This->stack[This->current - 1];
1129
1130     return D3D_OK;
1131 }
1132
1133 static HRESULT WINAPI ID3DXMatrixStackImpl_RotateAxis(ID3DXMatrixStack *iface, const D3DXVECTOR3 *pv, FLOAT angle)
1134 {
1135     D3DXMATRIX temp;
1136     struct ID3DXMatrixStackImpl *This = impl_from_ID3DXMatrixStack(iface);
1137
1138     TRACE("iface %p, pv %p, angle %f\n", iface, pv, angle);
1139
1140     D3DXMatrixRotationAxis(&temp, pv, angle);
1141     D3DXMatrixMultiply(&This->stack[This->current], &This->stack[This->current], &temp);
1142
1143     return D3D_OK;
1144 }
1145
1146 static HRESULT WINAPI ID3DXMatrixStackImpl_RotateAxisLocal(ID3DXMatrixStack *iface, const D3DXVECTOR3 *pv, FLOAT angle)
1147 {
1148     D3DXMATRIX temp;
1149     struct ID3DXMatrixStackImpl *This = impl_from_ID3DXMatrixStack(iface);
1150
1151     TRACE("iface %p, pv %p, angle %f\n", iface, pv, angle);
1152
1153     D3DXMatrixRotationAxis(&temp, pv, angle);
1154     D3DXMatrixMultiply(&This->stack[This->current], &temp, &This->stack[This->current]);
1155
1156     return D3D_OK;
1157 }
1158
1159 static HRESULT WINAPI ID3DXMatrixStackImpl_RotateYawPitchRoll(ID3DXMatrixStack *iface, FLOAT x, FLOAT y, FLOAT z)
1160 {
1161     D3DXMATRIX temp;
1162     struct ID3DXMatrixStackImpl *This = impl_from_ID3DXMatrixStack(iface);
1163
1164     TRACE("iface %p, x %f, y %f, z %f\n", iface, x, y, z);
1165
1166     D3DXMatrixRotationYawPitchRoll(&temp, x, y, z);
1167     D3DXMatrixMultiply(&This->stack[This->current], &This->stack[This->current], &temp);
1168
1169     return D3D_OK;
1170 }
1171
1172 static HRESULT WINAPI ID3DXMatrixStackImpl_RotateYawPitchRollLocal(ID3DXMatrixStack *iface, FLOAT x, FLOAT y, FLOAT z)
1173 {
1174     D3DXMATRIX temp;
1175     struct ID3DXMatrixStackImpl *This = impl_from_ID3DXMatrixStack(iface);
1176
1177     TRACE("iface %p, x %f, y %f, z %f\n", iface, x, y, z);
1178
1179     D3DXMatrixRotationYawPitchRoll(&temp, x, y, z);
1180     D3DXMatrixMultiply(&This->stack[This->current], &temp, &This->stack[This->current]);
1181
1182     return D3D_OK;
1183 }
1184
1185 static HRESULT WINAPI ID3DXMatrixStackImpl_Scale(ID3DXMatrixStack *iface, FLOAT x, FLOAT y, FLOAT z)
1186 {
1187     D3DXMATRIX temp;
1188     struct ID3DXMatrixStackImpl *This = impl_from_ID3DXMatrixStack(iface);
1189
1190     TRACE("iface %p,x %f, y %f, z %f\n", iface, x, y, z);
1191
1192     D3DXMatrixScaling(&temp, x, y, z);
1193     D3DXMatrixMultiply(&This->stack[This->current], &This->stack[This->current], &temp);
1194
1195     return D3D_OK;
1196 }
1197
1198 static HRESULT WINAPI ID3DXMatrixStackImpl_ScaleLocal(ID3DXMatrixStack *iface, FLOAT x, FLOAT y, FLOAT z)
1199 {
1200     D3DXMATRIX temp;
1201     struct ID3DXMatrixStackImpl *This = impl_from_ID3DXMatrixStack(iface);
1202
1203     TRACE("iface %p, x %f, y %f, z %f\n", iface, x, y, z);
1204
1205     D3DXMatrixScaling(&temp, x, y, z);
1206     D3DXMatrixMultiply(&This->stack[This->current], &temp, &This->stack[This->current]);
1207
1208     return D3D_OK;
1209 }
1210
1211 static HRESULT WINAPI ID3DXMatrixStackImpl_Translate(ID3DXMatrixStack *iface, FLOAT x, FLOAT y, FLOAT z)
1212 {
1213     D3DXMATRIX temp;
1214     struct ID3DXMatrixStackImpl *This = impl_from_ID3DXMatrixStack(iface);
1215
1216     TRACE("iface %p, x %f, y %f, z %f\n", iface, x, y, z);
1217
1218     D3DXMatrixTranslation(&temp, x, y, z);
1219     D3DXMatrixMultiply(&This->stack[This->current], &This->stack[This->current], &temp);
1220
1221     return D3D_OK;
1222 }
1223
1224 static HRESULT WINAPI ID3DXMatrixStackImpl_TranslateLocal(ID3DXMatrixStack *iface, FLOAT x, FLOAT y, FLOAT z)
1225 {
1226     D3DXMATRIX temp;
1227     struct ID3DXMatrixStackImpl *This = impl_from_ID3DXMatrixStack(iface);
1228
1229     TRACE("iface %p, x %f, y %f, z %f\n", iface, x, y, z);
1230
1231     D3DXMatrixTranslation(&temp, x, y, z);
1232     D3DXMatrixMultiply(&This->stack[This->current], &temp,&This->stack[This->current]);
1233
1234     return D3D_OK;
1235 }
1236
1237 static const ID3DXMatrixStackVtbl ID3DXMatrixStack_Vtbl =
1238 {
1239     ID3DXMatrixStackImpl_QueryInterface,
1240     ID3DXMatrixStackImpl_AddRef,
1241     ID3DXMatrixStackImpl_Release,
1242     ID3DXMatrixStackImpl_Pop,
1243     ID3DXMatrixStackImpl_Push,
1244     ID3DXMatrixStackImpl_LoadIdentity,
1245     ID3DXMatrixStackImpl_LoadMatrix,
1246     ID3DXMatrixStackImpl_MultMatrix,
1247     ID3DXMatrixStackImpl_MultMatrixLocal,
1248     ID3DXMatrixStackImpl_RotateAxis,
1249     ID3DXMatrixStackImpl_RotateAxisLocal,
1250     ID3DXMatrixStackImpl_RotateYawPitchRoll,
1251     ID3DXMatrixStackImpl_RotateYawPitchRollLocal,
1252     ID3DXMatrixStackImpl_Scale,
1253     ID3DXMatrixStackImpl_ScaleLocal,
1254     ID3DXMatrixStackImpl_Translate,
1255     ID3DXMatrixStackImpl_TranslateLocal,
1256     ID3DXMatrixStackImpl_GetTop
1257 };
1258
1259 /*_________________D3DXPLANE________________*/
1260
1261 D3DXPLANE* WINAPI D3DXPlaneFromPointNormal(D3DXPLANE *pout, const D3DXVECTOR3 *pvpoint, const D3DXVECTOR3 *pvnormal)
1262 {
1263     TRACE("pout %p, pvpoint %p, pvnormal %p\n", pout, pvpoint, pvnormal);
1264
1265     pout->a = pvnormal->x;
1266     pout->b = pvnormal->y;
1267     pout->c = pvnormal->z;
1268     pout->d = -D3DXVec3Dot(pvpoint, pvnormal);
1269     return pout;
1270 }
1271
1272 D3DXPLANE* WINAPI D3DXPlaneFromPoints(D3DXPLANE *pout, const D3DXVECTOR3 *pv1, const D3DXVECTOR3 *pv2, const D3DXVECTOR3 *pv3)
1273 {
1274     D3DXVECTOR3 edge1, edge2, normal, Nnormal;
1275
1276     TRACE("pout %p, pv1 %p, pv2 %p, pv3 %p\n", pout, pv1, pv2, pv3);
1277
1278     edge1.x = 0.0f; edge1.y = 0.0f; edge1.z = 0.0f;
1279     edge2.x = 0.0f; edge2.y = 0.0f; edge2.z = 0.0f;
1280     D3DXVec3Subtract(&edge1, pv2, pv1);
1281     D3DXVec3Subtract(&edge2, pv3, pv1);
1282     D3DXVec3Cross(&normal, &edge1, &edge2);
1283     D3DXVec3Normalize(&Nnormal, &normal);
1284     D3DXPlaneFromPointNormal(pout, pv1, &Nnormal);
1285     return pout;
1286 }
1287
1288 D3DXVECTOR3* WINAPI D3DXPlaneIntersectLine(D3DXVECTOR3 *pout, const D3DXPLANE *pp, const D3DXVECTOR3 *pv1, const D3DXVECTOR3 *pv2)
1289 {
1290     D3DXVECTOR3 direction, normal;
1291     FLOAT dot, temp;
1292
1293     TRACE("pout %p, pp %p, pv1 %p, pv2 %p\n", pout, pp, pv1, pv2);
1294
1295     normal.x = pp->a;
1296     normal.y = pp->b;
1297     normal.z = pp->c;
1298     direction.x = pv2->x - pv1->x;
1299     direction.y = pv2->y - pv1->y;
1300     direction.z = pv2->z - pv1->z;
1301     dot = D3DXVec3Dot(&normal, &direction);
1302     if ( !dot ) return NULL;
1303     temp = ( pp->d + D3DXVec3Dot(&normal, pv1) ) / dot;
1304     pout->x = pv1->x - temp * direction.x;
1305     pout->y = pv1->y - temp * direction.y;
1306     pout->z = pv1->z - temp * direction.z;
1307     return pout;
1308 }
1309
1310 D3DXPLANE * WINAPI D3DXPlaneNormalize(D3DXPLANE *out, const D3DXPLANE *p)
1311 {
1312     FLOAT norm;
1313
1314     TRACE("out %p, p %p\n", out, p);
1315
1316     norm = sqrtf(p->a * p->a + p->b * p->b + p->c * p->c);
1317     if (norm)
1318     {
1319         out->a = p->a / norm;
1320         out->b = p->b / norm;
1321         out->c = p->c / norm;
1322         out->d = p->d / norm;
1323     }
1324     else
1325     {
1326         out->a = 0.0f;
1327         out->b = 0.0f;
1328         out->c = 0.0f;
1329         out->d = 0.0f;
1330     }
1331
1332     return out;
1333 }
1334
1335 D3DXPLANE* WINAPI D3DXPlaneTransform(D3DXPLANE *pout, const D3DXPLANE *pplane, const D3DXMATRIX *pm)
1336 {
1337     const D3DXPLANE plane = *pplane;
1338
1339     TRACE("pout %p, pplane %p, pm %p\n", pout, pplane, pm);
1340
1341     pout->a = pm->u.m[0][0] * plane.a + pm->u.m[1][0] * plane.b + pm->u.m[2][0] * plane.c + pm->u.m[3][0] * plane.d;
1342     pout->b = pm->u.m[0][1] * plane.a + pm->u.m[1][1] * plane.b + pm->u.m[2][1] * plane.c + pm->u.m[3][1] * plane.d;
1343     pout->c = pm->u.m[0][2] * plane.a + pm->u.m[1][2] * plane.b + pm->u.m[2][2] * plane.c + pm->u.m[3][2] * plane.d;
1344     pout->d = pm->u.m[0][3] * plane.a + pm->u.m[1][3] * plane.b + pm->u.m[2][3] * plane.c + pm->u.m[3][3] * plane.d;
1345     return pout;
1346 }
1347
1348 D3DXPLANE* WINAPI D3DXPlaneTransformArray(D3DXPLANE* out, UINT outstride, const D3DXPLANE* in, UINT instride, const D3DXMATRIX* matrix, UINT elements)
1349 {
1350     UINT i;
1351
1352     TRACE("out %p, outstride %u, in %p, instride %u, matrix %p, elements %u\n", out, outstride, in, instride, matrix, elements);
1353
1354     for (i = 0; i < elements; ++i) {
1355         D3DXPlaneTransform(
1356             (D3DXPLANE*)((char*)out + outstride * i),
1357             (const D3DXPLANE*)((const char*)in + instride * i),
1358             matrix);
1359     }
1360     return out;
1361 }
1362
1363 /*_________________D3DXQUATERNION________________*/
1364
1365 D3DXQUATERNION* WINAPI D3DXQuaternionBaryCentric(D3DXQUATERNION *pout, const D3DXQUATERNION *pq1, const D3DXQUATERNION *pq2, const D3DXQUATERNION *pq3, FLOAT f, FLOAT g)
1366 {
1367     D3DXQUATERNION temp1, temp2;
1368
1369      TRACE("pout %p, pq1 %p, pq2 %p, pq3 %p, f %f, g %f\n", pout, pq1, pq2, pq3, f, g);
1370
1371     D3DXQuaternionSlerp(pout, D3DXQuaternionSlerp(&temp1, pq1, pq2, f + g), D3DXQuaternionSlerp(&temp2, pq1, pq3, f+g), g / (f + g));
1372     return pout;
1373 }
1374
1375 D3DXQUATERNION * WINAPI D3DXQuaternionExp(D3DXQUATERNION *out, const D3DXQUATERNION *q)
1376 {
1377     FLOAT norm;
1378
1379     TRACE("out %p, q %p\n", out, q);
1380
1381     norm = sqrtf(q->x * q->x + q->y * q->y + q->z * q->z);
1382     if (norm)
1383     {
1384         out->x = sinf(norm) * q->x / norm;
1385         out->y = sinf(norm) * q->y / norm;
1386         out->z = sinf(norm) * q->z / norm;
1387         out->w = cosf(norm);
1388     }
1389     else
1390     {
1391         out->x = 0.0f;
1392         out->y = 0.0f;
1393         out->z = 0.0f;
1394         out->w = 1.0f;
1395     }
1396
1397     return out;
1398 }
1399
1400 D3DXQUATERNION* WINAPI D3DXQuaternionInverse(D3DXQUATERNION *pout, const D3DXQUATERNION *pq)
1401 {
1402     D3DXQUATERNION out;
1403     FLOAT norm;
1404
1405     TRACE("pout %p, pq %p\n", pout, pq);
1406
1407     norm = D3DXQuaternionLengthSq(pq);
1408
1409     out.x = -pq->x / norm;
1410     out.y = -pq->y / norm;
1411     out.z = -pq->z / norm;
1412     out.w = pq->w / norm;
1413
1414     *pout =out;
1415     return pout;
1416 }
1417
1418 D3DXQUATERNION * WINAPI D3DXQuaternionLn(D3DXQUATERNION *out, const D3DXQUATERNION *q)
1419 {
1420     FLOAT t;
1421
1422     TRACE("out %p, q %p\n", out, q);
1423
1424     if ((q->w >= 1.0f) || (q->w == -1.0f))
1425         t = 1.0f;
1426     else
1427         t = acosf(q->w) / sqrtf(1.0f - q->w * q->w);
1428
1429     out->x = t * q->x;
1430     out->y = t * q->y;
1431     out->z = t * q->z;
1432     out->w = 0.0f;
1433
1434     return out;
1435 }
1436
1437 D3DXQUATERNION* WINAPI D3DXQuaternionMultiply(D3DXQUATERNION *pout, const D3DXQUATERNION *pq1, const D3DXQUATERNION *pq2)
1438 {
1439     D3DXQUATERNION out;
1440
1441     TRACE("pout %p, pq1 %p, pq2 %p\n", pout, pq1, pq2);
1442
1443     out.x = pq2->w * pq1->x + pq2->x * pq1->w + pq2->y * pq1->z - pq2->z * pq1->y;
1444     out.y = pq2->w * pq1->y - pq2->x * pq1->z + pq2->y * pq1->w + pq2->z * pq1->x;
1445     out.z = pq2->w * pq1->z + pq2->x * pq1->y - pq2->y * pq1->x + pq2->z * pq1->w;
1446     out.w = pq2->w * pq1->w - pq2->x * pq1->x - pq2->y * pq1->y - pq2->z * pq1->z;
1447     *pout = out;
1448     return pout;
1449 }
1450
1451 D3DXQUATERNION * WINAPI D3DXQuaternionNormalize(D3DXQUATERNION *out, const D3DXQUATERNION *q)
1452 {
1453     FLOAT norm;
1454
1455     TRACE("out %p, q %p\n", out, q);
1456
1457     norm = D3DXQuaternionLength(q);
1458
1459     out->x = q->x / norm;
1460     out->y = q->y / norm;
1461     out->z = q->z / norm;
1462     out->w = q->w / norm;
1463
1464     return out;
1465 }
1466
1467 D3DXQUATERNION * WINAPI D3DXQuaternionRotationAxis(D3DXQUATERNION *out, const D3DXVECTOR3 *v, FLOAT angle)
1468 {
1469     D3DXVECTOR3 temp;
1470
1471     TRACE("out %p, v %p, angle %f\n", out, v, angle);
1472
1473     D3DXVec3Normalize(&temp, v);
1474
1475     out->x = sinf(angle / 2.0f) * temp.x;
1476     out->y = sinf(angle / 2.0f) * temp.y;
1477     out->z = sinf(angle / 2.0f) * temp.z;
1478     out->w = cosf(angle / 2.0f);
1479
1480     return out;
1481 }
1482
1483 D3DXQUATERNION * WINAPI D3DXQuaternionRotationMatrix(D3DXQUATERNION *out, const D3DXMATRIX *m)
1484 {
1485     FLOAT s, trace;
1486
1487     TRACE("out %p, m %p\n", out, m);
1488
1489     trace = m->u.m[0][0] + m->u.m[1][1] + m->u.m[2][2] + 1.0f;
1490     if (trace > 1.0f)
1491     {
1492         s = 2.0f * sqrtf(trace);
1493         out->x = (m->u.m[1][2] - m->u.m[2][1]) / s;
1494         out->y = (m->u.m[2][0] - m->u.m[0][2]) / s;
1495         out->z = (m->u.m[0][1] - m->u.m[1][0]) / s;
1496         out->w = 0.25f * s;
1497     }
1498     else
1499     {
1500         int i, maxi = 0;
1501
1502         for (i = 1; i < 3; i++)
1503         {
1504             if (m->u.m[i][i] > m->u.m[maxi][maxi])
1505                 maxi = i;
1506         }
1507
1508         switch (maxi)
1509         {
1510             case 0:
1511                 s = 2.0f * sqrtf(1.0f + m->u.m[0][0] - m->u.m[1][1] - m->u.m[2][2]);
1512                 out->x = 0.25f * s;
1513                 out->y = (m->u.m[0][1] + m->u.m[1][0]) / s;
1514                 out->z = (m->u.m[0][2] + m->u.m[2][0]) / s;
1515                 out->w = (m->u.m[1][2] - m->u.m[2][1]) / s;
1516                 break;
1517
1518             case 1:
1519                 s = 2.0f * sqrtf(1.0f + m->u.m[1][1] - m->u.m[0][0] - m->u.m[2][2]);
1520                 out->x = (m->u.m[0][1] + m->u.m[1][0]) / s;
1521                 out->y = 0.25f * s;
1522                 out->z = (m->u.m[1][2] + m->u.m[2][1]) / s;
1523                 out->w = (m->u.m[2][0] - m->u.m[0][2]) / s;
1524                 break;
1525
1526             case 2:
1527                 s = 2.0f * sqrtf(1.0f + m->u.m[2][2] - m->u.m[0][0] - m->u.m[1][1]);
1528                 out->x = (m->u.m[0][2] + m->u.m[2][0]) / s;
1529                 out->y = (m->u.m[1][2] + m->u.m[2][1]) / s;
1530                 out->z = 0.25f * s;
1531                 out->w = (m->u.m[0][1] - m->u.m[1][0]) / s;
1532                 break;
1533         }
1534     }
1535
1536     return out;
1537 }
1538
1539 D3DXQUATERNION * WINAPI D3DXQuaternionRotationYawPitchRoll(D3DXQUATERNION *out, FLOAT yaw, FLOAT pitch, FLOAT roll)
1540 {
1541     FLOAT syaw, cyaw, spitch, cpitch, sroll, croll;
1542
1543     TRACE("out %p, yaw %f, pitch %f, roll %f\n", out, yaw, pitch, roll);
1544
1545     syaw = sinf(yaw / 2.0f);
1546     cyaw = cosf(yaw / 2.0f);
1547     spitch = sinf(pitch / 2.0f);
1548     cpitch = cosf(pitch / 2.0f);
1549     sroll = sinf(roll / 2.0f);
1550     croll = cosf(roll / 2.0f);
1551
1552     out->x = syaw * cpitch * sroll + cyaw * spitch * croll;
1553     out->y = syaw * cpitch * croll - cyaw * spitch * sroll;
1554     out->z = cyaw * cpitch * sroll - syaw * spitch * croll;
1555     out->w = cyaw * cpitch * croll + syaw * spitch * sroll;
1556
1557     return out;
1558 }
1559
1560 D3DXQUATERNION * WINAPI D3DXQuaternionSlerp(D3DXQUATERNION *out, const D3DXQUATERNION *q1,
1561         const D3DXQUATERNION *q2, FLOAT t)
1562 {
1563     FLOAT dot, temp;
1564
1565     TRACE("out %p, q1 %p, q2 %p, t %f\n", out, q1, q2, t);
1566
1567     temp = 1.0f - t;
1568     dot = D3DXQuaternionDot(q1, q2);
1569     if (dot < 0.0f)
1570     {
1571         t = -t;
1572         dot = -dot;
1573     }
1574
1575     if (1.0f - dot > 0.001f)
1576     {
1577         FLOAT theta = acosf(dot);
1578
1579         temp = sinf(theta * temp) / sinf(theta);
1580         t = sinf(theta * t) / sinf(theta);
1581     }
1582
1583     out->x = temp * q1->x + t * q2->x;
1584     out->y = temp * q1->y + t * q2->y;
1585     out->z = temp * q1->z + t * q2->z;
1586     out->w = temp * q1->w + t * q2->w;
1587
1588     return out;
1589 }
1590
1591 D3DXQUATERNION* WINAPI D3DXQuaternionSquad(D3DXQUATERNION *pout, const D3DXQUATERNION *pq1, const D3DXQUATERNION *pq2, const D3DXQUATERNION *pq3, const D3DXQUATERNION *pq4, FLOAT t)
1592 {
1593     D3DXQUATERNION temp1, temp2;
1594
1595     TRACE("pout %p, pq1 %p, pq2 %p, pq3 %p, pq4 %p, t %f\n", pout, pq1, pq2, pq3, pq4, t);
1596
1597     D3DXQuaternionSlerp(pout, D3DXQuaternionSlerp(&temp1, pq1, pq4, t), D3DXQuaternionSlerp(&temp2, pq2, pq3, t), 2.0f * t * (1.0f - t));
1598     return pout;
1599 }
1600
1601 static D3DXQUATERNION add_diff(const D3DXQUATERNION *q1, const D3DXQUATERNION *q2, const FLOAT add)
1602 {
1603     D3DXQUATERNION temp;
1604
1605     temp.x = q1->x + add * q2->x;
1606     temp.y = q1->y + add * q2->y;
1607     temp.z = q1->z + add * q2->z;
1608     temp.w = q1->w + add * q2->w;
1609
1610     return temp;
1611 }
1612
1613 void WINAPI D3DXQuaternionSquadSetup(D3DXQUATERNION *paout, D3DXQUATERNION *pbout, D3DXQUATERNION *pcout, const D3DXQUATERNION *pq0, const D3DXQUATERNION *pq1, const D3DXQUATERNION *pq2, const D3DXQUATERNION *pq3)
1614 {
1615     D3DXQUATERNION q, temp1, temp2, temp3, zero;
1616
1617     TRACE("paout %p, pbout %p, pcout %p, pq0 %p, pq1 %p, pq2 %p, pq3 %p\n", paout, pbout, pcout, pq0, pq1, pq2, pq3);
1618
1619     zero.x = 0.0f;
1620     zero.y = 0.0f;
1621     zero.z = 0.0f;
1622     zero.w = 0.0f;
1623
1624     if ( D3DXQuaternionDot(pq0, pq1) <  0.0f )
1625         temp2 = add_diff(&zero, pq0, -1.0f);
1626     else
1627         temp2 = *pq0;
1628
1629     if ( D3DXQuaternionDot(pq1, pq2) < 0.0f )
1630         *pcout = add_diff(&zero, pq2, -1.0f);
1631     else
1632         *pcout = *pq2;
1633
1634     if ( D3DXQuaternionDot(pcout, pq3) < 0.0f )
1635         temp3 = add_diff(&zero, pq3, -1.0f);
1636     else
1637         temp3 = *pq3;
1638
1639     D3DXQuaternionInverse(&temp1, pq1);
1640     D3DXQuaternionMultiply(&temp2, &temp1, &temp2);
1641     D3DXQuaternionLn(&temp2, &temp2);
1642     D3DXQuaternionMultiply(&q, &temp1, pcout);
1643     D3DXQuaternionLn(&q, &q);
1644     temp1 = add_diff(&temp2, &q, 1.0f);
1645     temp1.x *= -0.25f;
1646     temp1.y *= -0.25f;
1647     temp1.z *= -0.25f;
1648     temp1.w *= -0.25f;
1649     D3DXQuaternionExp(&temp1, &temp1);
1650     D3DXQuaternionMultiply(paout, pq1, &temp1);
1651
1652     D3DXQuaternionInverse(&temp1, pcout);
1653     D3DXQuaternionMultiply(&temp2, &temp1, pq1);
1654     D3DXQuaternionLn(&temp2, &temp2);
1655     D3DXQuaternionMultiply(&q, &temp1, &temp3);
1656     D3DXQuaternionLn(&q, &q);
1657     temp1 = add_diff(&temp2, &q, 1.0f);
1658     temp1.x *= -0.25f;
1659     temp1.y *= -0.25f;
1660     temp1.z *= -0.25f;
1661     temp1.w *= -0.25f;
1662     D3DXQuaternionExp(&temp1, &temp1);
1663     D3DXQuaternionMultiply(pbout, pcout, &temp1);
1664
1665     return;
1666 }
1667
1668 void WINAPI D3DXQuaternionToAxisAngle(const D3DXQUATERNION *pq, D3DXVECTOR3 *paxis, FLOAT *pangle)
1669 {
1670     TRACE("pq %p, paxis %p, pangle %p\n", pq, paxis, pangle);
1671
1672     paxis->x = pq->x;
1673     paxis->y = pq->y;
1674     paxis->z = pq->z;
1675     *pangle = 2.0f * acos(pq->w);
1676 }
1677
1678 /*_________________D3DXVec2_____________________*/
1679
1680 D3DXVECTOR2* WINAPI D3DXVec2BaryCentric(D3DXVECTOR2 *pout, const D3DXVECTOR2 *pv1, const D3DXVECTOR2 *pv2, const D3DXVECTOR2 *pv3, FLOAT f, FLOAT g)
1681 {
1682     TRACE("pout %p, pv1 %p, pv2 %p, pv3 %p, f %f, g %f\n", pout, pv1, pv2, pv3, f, g);
1683
1684     pout->x = (1.0f-f-g) * (pv1->x) + f * (pv2->x) + g * (pv3->x);
1685     pout->y = (1.0f-f-g) * (pv1->y) + f * (pv2->y) + g * (pv3->y);
1686     return pout;
1687 }
1688
1689 D3DXVECTOR2* WINAPI D3DXVec2CatmullRom(D3DXVECTOR2 *pout, const D3DXVECTOR2 *pv0, const D3DXVECTOR2 *pv1, const D3DXVECTOR2 *pv2, const D3DXVECTOR2 *pv3, FLOAT s)
1690 {
1691     TRACE("pout %p, pv0 %p, pv1 %p, pv2 %p, pv3 %p, s %f\n", pout, pv0, pv1, pv2, pv3, s);
1692
1693     pout->x = 0.5f * (2.0f * pv1->x + (pv2->x - pv0->x) *s + (2.0f *pv0->x - 5.0f * pv1->x + 4.0f * pv2->x - pv3->x) * s * s + (pv3->x -3.0f * pv2->x + 3.0f * pv1->x - pv0->x) * s * s * s);
1694     pout->y = 0.5f * (2.0f * pv1->y + (pv2->y - pv0->y) *s + (2.0f *pv0->y - 5.0f * pv1->y + 4.0f * pv2->y - pv3->y) * s * s + (pv3->y -3.0f * pv2->y + 3.0f * pv1->y - pv0->y) * s * s * s);
1695     return pout;
1696 }
1697
1698 D3DXVECTOR2* WINAPI D3DXVec2Hermite(D3DXVECTOR2 *pout, const D3DXVECTOR2 *pv1, const D3DXVECTOR2 *pt1, const D3DXVECTOR2 *pv2, const D3DXVECTOR2 *pt2, FLOAT s)
1699 {
1700     FLOAT h1, h2, h3, h4;
1701
1702     TRACE("pout %p, pv1 %p, pt1 %p, pv2 %p, pt2 %p, s %f\n", pout, pv1, pt1, pv2, pt2, s);
1703
1704     h1 = 2.0f * s * s * s - 3.0f * s * s + 1.0f;
1705     h2 = s * s * s - 2.0f * s * s + s;
1706     h3 = -2.0f * s * s * s + 3.0f * s * s;
1707     h4 = s * s * s - s * s;
1708
1709     pout->x = h1 * (pv1->x) + h2 * (pt1->x) + h3 * (pv2->x) + h4 * (pt2->x);
1710     pout->y = h1 * (pv1->y) + h2 * (pt1->y) + h3 * (pv2->y) + h4 * (pt2->y);
1711     return pout;
1712 }
1713
1714 D3DXVECTOR2* WINAPI D3DXVec2Normalize(D3DXVECTOR2 *pout, const D3DXVECTOR2 *pv)
1715 {
1716     FLOAT norm;
1717
1718     TRACE("pout %p, pv %p\n", pout, pv);
1719
1720     norm = D3DXVec2Length(pv);
1721     if ( !norm )
1722     {
1723         pout->x = 0.0f;
1724         pout->y = 0.0f;
1725     }
1726     else
1727     {
1728         pout->x = pv->x / norm;
1729         pout->y = pv->y / norm;
1730     }
1731
1732     return pout;
1733 }
1734
1735 D3DXVECTOR4* WINAPI D3DXVec2Transform(D3DXVECTOR4 *pout, const D3DXVECTOR2 *pv, const D3DXMATRIX *pm)
1736 {
1737     TRACE("pout %p, pv %p, pm %p\n", pout, pv, pm);
1738
1739     pout->x = pm->u.m[0][0] * pv->x + pm->u.m[1][0] * pv->y  + pm->u.m[3][0];
1740     pout->y = pm->u.m[0][1] * pv->x + pm->u.m[1][1] * pv->y  + pm->u.m[3][1];
1741     pout->z = pm->u.m[0][2] * pv->x + pm->u.m[1][2] * pv->y  + pm->u.m[3][2];
1742     pout->w = pm->u.m[0][3] * pv->x + pm->u.m[1][3] * pv->y  + pm->u.m[3][3];
1743     return pout;
1744 }
1745
1746 D3DXVECTOR4* WINAPI D3DXVec2TransformArray(D3DXVECTOR4* out, UINT outstride, const D3DXVECTOR2* in, UINT instride, const D3DXMATRIX* matrix, UINT elements)
1747 {
1748     UINT i;
1749
1750     TRACE("out %p, outstride %u, in %p, instride %u, matrix %p, elements %u\n", out, outstride, in, instride, matrix, elements);
1751
1752     for (i = 0; i < elements; ++i) {
1753         D3DXVec2Transform(
1754             (D3DXVECTOR4*)((char*)out + outstride * i),
1755             (const D3DXVECTOR2*)((const char*)in + instride * i),
1756             matrix);
1757     }
1758     return out;
1759 }
1760
1761 D3DXVECTOR2* WINAPI D3DXVec2TransformCoord(D3DXVECTOR2 *pout, const D3DXVECTOR2 *pv, const D3DXMATRIX *pm)
1762 {
1763     D3DXVECTOR2 v;
1764     FLOAT norm;
1765
1766     TRACE("pout %p, pv %p, pm %p\n", pout, pv, pm);
1767
1768     v = *pv;
1769     norm = pm->u.m[0][3] * pv->x + pm->u.m[1][3] * pv->y + pm->u.m[3][3];
1770
1771     pout->x = (pm->u.m[0][0] * v.x + pm->u.m[1][0] * v.y + pm->u.m[3][0]) / norm;
1772     pout->y = (pm->u.m[0][1] * v.x + pm->u.m[1][1] * v.y + pm->u.m[3][1]) / norm;
1773
1774     return pout;
1775 }
1776
1777 D3DXVECTOR2* WINAPI D3DXVec2TransformCoordArray(D3DXVECTOR2* out, UINT outstride, const D3DXVECTOR2* in, UINT instride, const D3DXMATRIX* matrix, UINT elements)
1778 {
1779     UINT i;
1780
1781     TRACE("out %p, outstride %u, in %p, instride %u, matrix %p, elements %u\n", out, outstride, in, instride, matrix, elements);
1782
1783     for (i = 0; i < elements; ++i) {
1784         D3DXVec2TransformCoord(
1785             (D3DXVECTOR2*)((char*)out + outstride * i),
1786             (const D3DXVECTOR2*)((const char*)in + instride * i),
1787             matrix);
1788     }
1789     return out;
1790 }
1791
1792 D3DXVECTOR2* WINAPI D3DXVec2TransformNormal(D3DXVECTOR2 *pout, const D3DXVECTOR2 *pv, const D3DXMATRIX *pm)
1793 {
1794     const D3DXVECTOR2 v = *pv;
1795
1796     TRACE("pout %p, pv %p, pm %p\n", pout, pv, pm);
1797
1798     pout->x = pm->u.m[0][0] * v.x + pm->u.m[1][0] * v.y;
1799     pout->y = pm->u.m[0][1] * v.x + pm->u.m[1][1] * v.y;
1800     return pout;
1801 }
1802
1803 D3DXVECTOR2* WINAPI D3DXVec2TransformNormalArray(D3DXVECTOR2* out, UINT outstride, const D3DXVECTOR2 *in, UINT instride, const D3DXMATRIX *matrix, UINT elements)
1804 {
1805     UINT i;
1806
1807     TRACE("out %p, outstride %u, in %p, instride %u, matrix %p, elements %u\n", out, outstride, in, instride, matrix, elements);
1808
1809     for (i = 0; i < elements; ++i) {
1810         D3DXVec2TransformNormal(
1811             (D3DXVECTOR2*)((char*)out + outstride * i),
1812             (const D3DXVECTOR2*)((const char*)in + instride * i),
1813             matrix);
1814     }
1815     return out;
1816 }
1817
1818 /*_________________D3DXVec3_____________________*/
1819
1820 D3DXVECTOR3* WINAPI D3DXVec3BaryCentric(D3DXVECTOR3 *pout, const D3DXVECTOR3 *pv1, const D3DXVECTOR3 *pv2, const D3DXVECTOR3 *pv3, FLOAT f, FLOAT g)
1821 {
1822     TRACE("pout %p, pv1 %p, pv2 %p, pv3 %p, f %f, g %f\n", pout, pv1, pv2, pv3, f, g);
1823
1824     pout->x = (1.0f-f-g) * (pv1->x) + f * (pv2->x) + g * (pv3->x);
1825     pout->y = (1.0f-f-g) * (pv1->y) + f * (pv2->y) + g * (pv3->y);
1826     pout->z = (1.0f-f-g) * (pv1->z) + f * (pv2->z) + g * (pv3->z);
1827     return pout;
1828 }
1829
1830 D3DXVECTOR3* WINAPI D3DXVec3CatmullRom( D3DXVECTOR3 *pout, const D3DXVECTOR3 *pv0, const D3DXVECTOR3 *pv1, const D3DXVECTOR3 *pv2, const D3DXVECTOR3 *pv3, FLOAT s)
1831 {
1832     TRACE("pout %p, pv0 %p, pv1 %p, pv2 %p, pv3 %p, s %f\n", pout, pv0, pv1, pv2, pv3, s);
1833
1834     pout->x = 0.5f * (2.0f * pv1->x + (pv2->x - pv0->x) *s + (2.0f *pv0->x - 5.0f * pv1->x + 4.0f * pv2->x - pv3->x) * s * s + (pv3->x -3.0f * pv2->x + 3.0f * pv1->x - pv0->x) * s * s * s);
1835     pout->y = 0.5f * (2.0f * pv1->y + (pv2->y - pv0->y) *s + (2.0f *pv0->y - 5.0f * pv1->y + 4.0f * pv2->y - pv3->y) * s * s + (pv3->y -3.0f * pv2->y + 3.0f * pv1->y - pv0->y) * s * s * s);
1836     pout->z = 0.5f * (2.0f * pv1->z + (pv2->z - pv0->z) *s + (2.0f *pv0->z - 5.0f * pv1->z + 4.0f * pv2->z - pv3->z) * s * s + (pv3->z -3.0f * pv2->z + 3.0f * pv1->z - pv0->z) * s * s * s);
1837     return pout;
1838 }
1839
1840 D3DXVECTOR3* WINAPI D3DXVec3Hermite(D3DXVECTOR3 *pout, const D3DXVECTOR3 *pv1, const D3DXVECTOR3 *pt1, const D3DXVECTOR3 *pv2, const D3DXVECTOR3 *pt2, FLOAT s)
1841 {
1842     FLOAT h1, h2, h3, h4;
1843
1844     TRACE("pout %p, pv1 %p, pt1 %p, pv2 %p, pt2 %p, s %f\n", pout, pv1, pt1, pv2, pt2, s);
1845
1846     h1 = 2.0f * s * s * s - 3.0f * s * s + 1.0f;
1847     h2 = s * s * s - 2.0f * s * s + s;
1848     h3 = -2.0f * s * s * s + 3.0f * s * s;
1849     h4 = s * s * s - s * s;
1850
1851     pout->x = h1 * (pv1->x) + h2 * (pt1->x) + h3 * (pv2->x) + h4 * (pt2->x);
1852     pout->y = h1 * (pv1->y) + h2 * (pt1->y) + h3 * (pv2->y) + h4 * (pt2->y);
1853     pout->z = h1 * (pv1->z) + h2 * (pt1->z) + h3 * (pv2->z) + h4 * (pt2->z);
1854     return pout;
1855 }
1856
1857 D3DXVECTOR3* WINAPI D3DXVec3Normalize(D3DXVECTOR3 *pout, const D3DXVECTOR3 *pv)
1858 {
1859     FLOAT norm;
1860
1861     TRACE("pout %p, pv %p\n", pout, pv);
1862
1863     norm = D3DXVec3Length(pv);
1864     if ( !norm )
1865     {
1866         pout->x = 0.0f;
1867         pout->y = 0.0f;
1868         pout->z = 0.0f;
1869     }
1870     else
1871     {
1872         pout->x = pv->x / norm;
1873         pout->y = pv->y / norm;
1874         pout->z = pv->z / norm;
1875     }
1876
1877     return pout;
1878 }
1879
1880 D3DXVECTOR3* WINAPI D3DXVec3Project(D3DXVECTOR3 *pout, const D3DXVECTOR3 *pv, const D3DVIEWPORT9 *pviewport, const D3DXMATRIX *pprojection, const D3DXMATRIX *pview, const D3DXMATRIX *pworld)
1881 {
1882     D3DXMATRIX m;
1883
1884     TRACE("pout %p, pv %p, pviewport %p, pprojection %p, pview %p, pworld %p\n", pout, pv, pviewport, pprojection, pview, pworld);
1885
1886     D3DXMatrixIdentity(&m);
1887     if (pworld) D3DXMatrixMultiply(&m, &m, pworld);
1888     if (pview) D3DXMatrixMultiply(&m, &m, pview);
1889     if (pprojection) D3DXMatrixMultiply(&m, &m, pprojection);
1890
1891     D3DXVec3TransformCoord(pout, pv, &m);
1892
1893     if (pviewport)
1894     {
1895         pout->x = pviewport->X +  ( 1.0f + pout->x ) * pviewport->Width / 2.0f;
1896         pout->y = pviewport->Y +  ( 1.0f - pout->y ) * pviewport->Height / 2.0f;
1897         pout->z = pviewport->MinZ + pout->z * ( pviewport->MaxZ - pviewport->MinZ );
1898     }
1899     return pout;
1900 }
1901
1902 D3DXVECTOR3* WINAPI D3DXVec3ProjectArray(D3DXVECTOR3* out, UINT outstride, const D3DXVECTOR3* in, UINT instride, const D3DVIEWPORT9* viewport, const D3DXMATRIX* projection, const D3DXMATRIX* view, const D3DXMATRIX* world, UINT elements)
1903 {
1904     UINT i;
1905
1906     TRACE("out %p, outstride %u, in %p, instride %u, viewport %p, projection %p, view %p, world %p, elements %u\n",
1907         out, outstride, in, instride, viewport, projection, view, world, elements);
1908
1909     for (i = 0; i < elements; ++i) {
1910         D3DXVec3Project(
1911             (D3DXVECTOR3*)((char*)out + outstride * i),
1912             (const D3DXVECTOR3*)((const char*)in + instride * i),
1913             viewport, projection, view, world);
1914     }
1915     return out;
1916 }
1917
1918 D3DXVECTOR4* WINAPI D3DXVec3Transform(D3DXVECTOR4 *pout, const D3DXVECTOR3 *pv, const D3DXMATRIX *pm)
1919 {
1920     TRACE("pout %p, pv %p, pm %p\n", pout, pv, pm);
1921
1922     pout->x = pm->u.m[0][0] * pv->x + pm->u.m[1][0] * pv->y + pm->u.m[2][0] * pv->z + pm->u.m[3][0];
1923     pout->y = pm->u.m[0][1] * pv->x + pm->u.m[1][1] * pv->y + pm->u.m[2][1] * pv->z + pm->u.m[3][1];
1924     pout->z = pm->u.m[0][2] * pv->x + pm->u.m[1][2] * pv->y + pm->u.m[2][2] * pv->z + pm->u.m[3][2];
1925     pout->w = pm->u.m[0][3] * pv->x + pm->u.m[1][3] * pv->y + pm->u.m[2][3] * pv->z + pm->u.m[3][3];
1926     return pout;
1927 }
1928
1929 D3DXVECTOR4* WINAPI D3DXVec3TransformArray(D3DXVECTOR4* out, UINT outstride, const D3DXVECTOR3* in, UINT instride, const D3DXMATRIX* matrix, UINT elements)
1930 {
1931     UINT i;
1932
1933     TRACE("out %p, outstride %u, in %p, instride %u, matrix %p, elements %u\n", out, outstride, in, instride, matrix, elements);
1934
1935     for (i = 0; i < elements; ++i) {
1936         D3DXVec3Transform(
1937             (D3DXVECTOR4*)((char*)out + outstride * i),
1938             (const D3DXVECTOR3*)((const char*)in + instride * i),
1939             matrix);
1940     }
1941     return out;
1942 }
1943
1944 D3DXVECTOR3* WINAPI D3DXVec3TransformCoord(D3DXVECTOR3 *pout, const D3DXVECTOR3 *pv, const D3DXMATRIX *pm)
1945 {
1946     D3DXVECTOR3 out;
1947     FLOAT norm;
1948
1949     TRACE("pout %p, pv %p, pm %p\n", pout, pv, pm);
1950
1951     norm = pm->u.m[0][3] * pv->x + pm->u.m[1][3] * pv->y + pm->u.m[2][3] *pv->z + pm->u.m[3][3];
1952
1953     out.x = (pm->u.m[0][0] * pv->x + pm->u.m[1][0] * pv->y + pm->u.m[2][0] * pv->z + pm->u.m[3][0]) / norm;
1954     out.y = (pm->u.m[0][1] * pv->x + pm->u.m[1][1] * pv->y + pm->u.m[2][1] * pv->z + pm->u.m[3][1]) / norm;
1955     out.z = (pm->u.m[0][2] * pv->x + pm->u.m[1][2] * pv->y + pm->u.m[2][2] * pv->z + pm->u.m[3][2]) / norm;
1956
1957     *pout = out;
1958
1959     return pout;
1960 }
1961
1962 D3DXVECTOR3* WINAPI D3DXVec3TransformCoordArray(D3DXVECTOR3* out, UINT outstride, const D3DXVECTOR3* in, UINT instride, const D3DXMATRIX* matrix, UINT elements)
1963 {
1964     UINT i;
1965
1966     TRACE("out %p, outstride %u, in %p, instride %u, matrix %p, elements %u\n", out, outstride, in, instride, matrix, elements);
1967
1968     for (i = 0; i < elements; ++i) {
1969         D3DXVec3TransformCoord(
1970             (D3DXVECTOR3*)((char*)out + outstride * i),
1971             (const D3DXVECTOR3*)((const char*)in + instride * i),
1972             matrix);
1973     }
1974     return out;
1975 }
1976
1977 D3DXVECTOR3* WINAPI D3DXVec3TransformNormal(D3DXVECTOR3 *pout, const D3DXVECTOR3 *pv, const D3DXMATRIX *pm)
1978 {
1979     const D3DXVECTOR3 v = *pv;
1980
1981     TRACE("pout %p, pv %p, pm %p\n", pout, pv, pm);
1982
1983     pout->x = pm->u.m[0][0] * v.x + pm->u.m[1][0] * v.y + pm->u.m[2][0] * v.z;
1984     pout->y = pm->u.m[0][1] * v.x + pm->u.m[1][1] * v.y + pm->u.m[2][1] * v.z;
1985     pout->z = pm->u.m[0][2] * v.x + pm->u.m[1][2] * v.y + pm->u.m[2][2] * v.z;
1986     return pout;
1987
1988 }
1989
1990 D3DXVECTOR3* WINAPI D3DXVec3TransformNormalArray(D3DXVECTOR3* out, UINT outstride, const D3DXVECTOR3* in, UINT instride, const D3DXMATRIX* matrix, UINT elements)
1991 {
1992     UINT i;
1993
1994     TRACE("out %p, outstride %u, in %p, instride %u, matrix %p, elements %u\n", out, outstride, in, instride, matrix, elements);
1995
1996     for (i = 0; i < elements; ++i) {
1997         D3DXVec3TransformNormal(
1998             (D3DXVECTOR3*)((char*)out + outstride * i),
1999             (const D3DXVECTOR3*)((const char*)in + instride * i),
2000             matrix);
2001     }
2002     return out;
2003 }
2004
2005 D3DXVECTOR3* WINAPI D3DXVec3Unproject(D3DXVECTOR3 *pout, const D3DXVECTOR3 *pv, const D3DVIEWPORT9 *pviewport, const D3DXMATRIX *pprojection, const D3DXMATRIX *pview, const D3DXMATRIX *pworld)
2006 {
2007     D3DXMATRIX m;
2008
2009     TRACE("pout %p, pv %p, pviewport %p, pprojection %p, pview %p, pworlds %p\n", pout, pv, pviewport, pprojection, pview, pworld);
2010
2011     D3DXMatrixIdentity(&m);
2012     if (pworld) D3DXMatrixMultiply(&m, &m, pworld);
2013     if (pview) D3DXMatrixMultiply(&m, &m, pview);
2014     if (pprojection) D3DXMatrixMultiply(&m, &m, pprojection);
2015     D3DXMatrixInverse(&m, NULL, &m);
2016
2017     *pout = *pv;
2018     if (pviewport)
2019     {
2020         pout->x = 2.0f * ( pout->x - pviewport->X ) / pviewport->Width - 1.0f;
2021         pout->y = 1.0f - 2.0f * ( pout->y - pviewport->Y ) / pviewport->Height;
2022         pout->z = ( pout->z - pviewport->MinZ) / ( pviewport->MaxZ - pviewport->MinZ );
2023     }
2024     D3DXVec3TransformCoord(pout, pout, &m);
2025     return pout;
2026 }
2027
2028 D3DXVECTOR3* WINAPI D3DXVec3UnprojectArray(D3DXVECTOR3* out, UINT outstride, const D3DXVECTOR3* in, UINT instride, const D3DVIEWPORT9* viewport, const D3DXMATRIX* projection, const D3DXMATRIX* view, const D3DXMATRIX* world, UINT elements)
2029 {
2030     UINT i;
2031
2032     TRACE("out %p, outstride %u, in %p, instride %u, viewport %p, projection %p, view %p, world %p, elements %u\n",
2033         out, outstride, in, instride, viewport, projection, view, world, elements);
2034
2035     for (i = 0; i < elements; ++i) {
2036         D3DXVec3Unproject(
2037             (D3DXVECTOR3*)((char*)out + outstride * i),
2038             (const D3DXVECTOR3*)((const char*)in + instride * i),
2039             viewport, projection, view, world);
2040     }
2041     return out;
2042 }
2043
2044 /*_________________D3DXVec4_____________________*/
2045
2046 D3DXVECTOR4* WINAPI D3DXVec4BaryCentric(D3DXVECTOR4 *pout, const D3DXVECTOR4 *pv1, const D3DXVECTOR4 *pv2, const D3DXVECTOR4 *pv3, FLOAT f, FLOAT g)
2047 {
2048     TRACE("pout %p, pv1 %p, pv2 %p, pv3 %p, f %f, g %f\n", pout, pv1, pv2, pv3, f, g);
2049
2050     pout->x = (1.0f-f-g) * (pv1->x) + f * (pv2->x) + g * (pv3->x);
2051     pout->y = (1.0f-f-g) * (pv1->y) + f * (pv2->y) + g * (pv3->y);
2052     pout->z = (1.0f-f-g) * (pv1->z) + f * (pv2->z) + g * (pv3->z);
2053     pout->w = (1.0f-f-g) * (pv1->w) + f * (pv2->w) + g * (pv3->w);
2054     return pout;
2055 }
2056
2057 D3DXVECTOR4* WINAPI D3DXVec4CatmullRom(D3DXVECTOR4 *pout, const D3DXVECTOR4 *pv0, const D3DXVECTOR4 *pv1, const D3DXVECTOR4 *pv2, const D3DXVECTOR4 *pv3, FLOAT s)
2058 {
2059     TRACE("pout %p, pv0 %p, pv1 %p, pv2 %p, pv3 %p, s %f\n", pout, pv0, pv1, pv2, pv3, s);
2060
2061     pout->x = 0.5f * (2.0f * pv1->x + (pv2->x - pv0->x) *s + (2.0f *pv0->x - 5.0f * pv1->x + 4.0f * pv2->x - pv3->x) * s * s + (pv3->x -3.0f * pv2->x + 3.0f * pv1->x - pv0->x) * s * s * s);
2062     pout->y = 0.5f * (2.0f * pv1->y + (pv2->y - pv0->y) *s + (2.0f *pv0->y - 5.0f * pv1->y + 4.0f * pv2->y - pv3->y) * s * s + (pv3->y -3.0f * pv2->y + 3.0f * pv1->y - pv0->y) * s * s * s);
2063     pout->z = 0.5f * (2.0f * pv1->z + (pv2->z - pv0->z) *s + (2.0f *pv0->z - 5.0f * pv1->z + 4.0f * pv2->z - pv3->z) * s * s + (pv3->z -3.0f * pv2->z + 3.0f * pv1->z - pv0->z) * s * s * s);
2064     pout->w = 0.5f * (2.0f * pv1->w + (pv2->w - pv0->w) *s + (2.0f *pv0->w - 5.0f * pv1->w + 4.0f * pv2->w - pv3->w) * s * s + (pv3->w -3.0f * pv2->w + 3.0f * pv1->w - pv0->w) * s * s * s);
2065     return pout;
2066 }
2067
2068 D3DXVECTOR4* WINAPI D3DXVec4Cross(D3DXVECTOR4 *pout, const D3DXVECTOR4 *pv1, const D3DXVECTOR4 *pv2, const D3DXVECTOR4 *pv3)
2069 {
2070     D3DXVECTOR4 out;
2071
2072     TRACE("pout %p, pv1 %p, pv2 %p, pv3 %p\n", pout, pv1, pv2, pv3);
2073
2074     out.x = pv1->y * (pv2->z * pv3->w - pv3->z * pv2->w) - pv1->z * (pv2->y * pv3->w - pv3->y * pv2->w) + pv1->w * (pv2->y * pv3->z - pv2->z *pv3->y);
2075     out.y = -(pv1->x * (pv2->z * pv3->w - pv3->z * pv2->w) - pv1->z * (pv2->x * pv3->w - pv3->x * pv2->w) + pv1->w * (pv2->x * pv3->z - pv3->x * pv2->z));
2076     out.z = pv1->x * (pv2->y * pv3->w - pv3->y * pv2->w) - pv1->y * (pv2->x *pv3->w - pv3->x * pv2->w) + pv1->w * (pv2->x * pv3->y - pv3->x * pv2->y);
2077     out.w = -(pv1->x * (pv2->y * pv3->z - pv3->y * pv2->z) - pv1->y * (pv2->x * pv3->z - pv3->x *pv2->z) + pv1->z * (pv2->x * pv3->y - pv3->x * pv2->y));
2078     *pout = out;
2079     return pout;
2080 }
2081
2082 D3DXVECTOR4* WINAPI D3DXVec4Hermite(D3DXVECTOR4 *pout, const D3DXVECTOR4 *pv1, const D3DXVECTOR4 *pt1, const D3DXVECTOR4 *pv2, const D3DXVECTOR4 *pt2, FLOAT s)
2083 {
2084     FLOAT h1, h2, h3, h4;
2085
2086     TRACE("pout %p, pv1 %p, pt1 %p, pv2 %p, pt2 %p, s %f\n", pout, pv1, pt1, pv2, pt2, s);
2087
2088     h1 = 2.0f * s * s * s - 3.0f * s * s + 1.0f;
2089     h2 = s * s * s - 2.0f * s * s + s;
2090     h3 = -2.0f * s * s * s + 3.0f * s * s;
2091     h4 = s * s * s - s * s;
2092
2093     pout->x = h1 * (pv1->x) + h2 * (pt1->x) + h3 * (pv2->x) + h4 * (pt2->x);
2094     pout->y = h1 * (pv1->y) + h2 * (pt1->y) + h3 * (pv2->y) + h4 * (pt2->y);
2095     pout->z = h1 * (pv1->z) + h2 * (pt1->z) + h3 * (pv2->z) + h4 * (pt2->z);
2096     pout->w = h1 * (pv1->w) + h2 * (pt1->w) + h3 * (pv2->w) + h4 * (pt2->w);
2097     return pout;
2098 }
2099
2100 D3DXVECTOR4* WINAPI D3DXVec4Normalize(D3DXVECTOR4 *pout, const D3DXVECTOR4 *pv)
2101 {
2102     FLOAT norm;
2103
2104     TRACE("pout %p, pv %p\n", pout, pv);
2105
2106     norm = D3DXVec4Length(pv);
2107
2108     pout->x = pv->x / norm;
2109     pout->y = pv->y / norm;
2110     pout->z = pv->z / norm;
2111     pout->w = pv->w / norm;
2112
2113     return pout;
2114 }
2115
2116 D3DXVECTOR4* WINAPI D3DXVec4Transform(D3DXVECTOR4 *pout, const D3DXVECTOR4 *pv, const D3DXMATRIX *pm)
2117 {
2118     D3DXVECTOR4 out;
2119
2120     TRACE("pout %p, pv %p, pm %p\n", pout, pv, pm);
2121
2122     out.x = pm->u.m[0][0] * pv->x + pm->u.m[1][0] * pv->y + pm->u.m[2][0] * pv->z + pm->u.m[3][0] * pv->w;
2123     out.y = pm->u.m[0][1] * pv->x + pm->u.m[1][1] * pv->y + pm->u.m[2][1] * pv->z + pm->u.m[3][1] * pv->w;
2124     out.z = pm->u.m[0][2] * pv->x + pm->u.m[1][2] * pv->y + pm->u.m[2][2] * pv->z + pm->u.m[3][2] * pv->w;
2125     out.w = pm->u.m[0][3] * pv->x + pm->u.m[1][3] * pv->y + pm->u.m[2][3] * pv->z + pm->u.m[3][3] * pv->w;
2126     *pout = out;
2127     return pout;
2128 }
2129
2130 D3DXVECTOR4* WINAPI D3DXVec4TransformArray(D3DXVECTOR4* out, UINT outstride, const D3DXVECTOR4* in, UINT instride, const D3DXMATRIX* matrix, UINT elements)
2131 {
2132     UINT i;
2133
2134     TRACE("out %p, outstride %u, in %p, instride %u, matrix %p, elements %u\n", out, outstride, in, instride, matrix, elements);
2135
2136     for (i = 0; i < elements; ++i) {
2137         D3DXVec4Transform(
2138             (D3DXVECTOR4*)((char*)out + outstride * i),
2139             (const D3DXVECTOR4*)((const char*)in + instride * i),
2140             matrix);
2141     }
2142     return out;
2143 }
2144
2145 unsigned short float_32_to_16(const float in)
2146 {
2147     int exp = 0, origexp;
2148     float tmp = fabs(in);
2149     int sign = (copysignf(1, in) < 0);
2150     unsigned int mantissa;
2151     unsigned short ret;
2152
2153     /* Deal with special numbers */
2154     if (isinf(in)) return (sign ? 0xffff : 0x7fff);
2155     if (isnan(in)) return (sign ? 0xffff : 0x7fff);
2156     if (in == 0.0f) return (sign ? 0x8000 : 0x0000);
2157
2158     if (tmp < powf(2, 10))
2159     {
2160         do
2161         {
2162             tmp *= 2.0f;
2163             exp--;
2164         } while (tmp < powf(2, 10));
2165     }
2166     else if (tmp >= powf(2, 11))
2167     {
2168         do
2169         {
2170             tmp /= 2.0f;
2171             exp++;
2172         } while (tmp >= powf(2, 11));
2173     }
2174
2175     exp += 10;  /* Normalize the mantissa */
2176     exp += 15;  /* Exponent is encoded with excess 15 */
2177
2178     origexp = exp;
2179
2180     mantissa = (unsigned int) tmp;
2181     if ((tmp - mantissa == 0.5f && mantissa % 2 == 1) || /* round half to even */
2182         (tmp - mantissa > 0.5f))
2183     {
2184         mantissa++; /* round to nearest, away from zero */
2185     }
2186     if (mantissa == 2048)
2187     {
2188         mantissa = 1024;
2189         exp++;
2190     }
2191
2192     if (exp > 31)
2193     {
2194         /* too big */
2195         ret = 0x7fff; /* INF */
2196     }
2197     else if (exp <= 0)
2198     {
2199         unsigned int rounding = 0;
2200
2201         /* Denormalized half float */
2202
2203         /* return 0x0000 (=0.0) for numbers too small to represent in half floats */
2204         if (exp < -11)
2205             return (sign ? 0x8000 : 0x0000);
2206
2207         exp = origexp;
2208
2209         /* the 13 extra bits from single precision are used for rounding */
2210         mantissa = (unsigned int)(tmp * powf(2, 13));
2211         mantissa >>= 1 - exp; /* denormalize */
2212
2213         mantissa -= ~(mantissa >> 13) & 1; /* round half to even */
2214         /* remove 13 least significant bits to get half float precision */
2215         mantissa >>= 12;
2216         rounding = mantissa & 1;
2217         mantissa >>= 1;
2218
2219         ret = mantissa + rounding;
2220     }
2221     else
2222     {
2223         ret = (exp << 10) | (mantissa & 0x3ff);
2224     }
2225
2226     ret |= ((sign ? 1 : 0) << 15); /* Add the sign */
2227     return ret;
2228 }
2229
2230 D3DXFLOAT16 *WINAPI D3DXFloat32To16Array(D3DXFLOAT16 *pout, const FLOAT *pin, UINT n)
2231 {
2232     unsigned int i;
2233
2234     TRACE("pout %p, pin %p, n %u\n", pout, pin, n);
2235
2236     for (i = 0; i < n; ++i)
2237     {
2238         pout[i].value = float_32_to_16(pin[i]);
2239     }
2240
2241     return pout;
2242 }
2243
2244 /* Native d3dx9's D3DXFloat16to32Array lacks support for NaN and Inf. Specifically, e = 16 is treated as a
2245  * regular number - e.g., 0x7fff is converted to 131008.0 and 0xffff to -131008.0. */
2246 static inline float float_16_to_32(const unsigned short in)
2247 {
2248     const unsigned short s = (in & 0x8000);
2249     const unsigned short e = (in & 0x7C00) >> 10;
2250     const unsigned short m = in & 0x3FF;
2251     const float sgn = (s ? -1.0f : 1.0f);
2252
2253     if (e == 0)
2254     {
2255         if (m == 0) return sgn * 0.0f; /* +0.0 or -0.0 */
2256         else return sgn * powf(2, -14.0f) * (m / 1024.0f);
2257     }
2258     else
2259     {
2260         return sgn * powf(2, e - 15.0f) * (1.0f + (m / 1024.0f));
2261     }
2262 }
2263
2264 FLOAT *WINAPI D3DXFloat16To32Array(FLOAT *pout, const D3DXFLOAT16 *pin, UINT n)
2265 {
2266     unsigned int i;
2267
2268     TRACE("pout %p, pin %p, n %u\n", pout, pin, n);
2269
2270     for (i = 0; i < n; ++i)
2271     {
2272         pout[i] = float_16_to_32(pin[i].value);
2273     }
2274
2275     return pout;
2276 }
2277
2278 /*_________________D3DXSH________________*/
2279
2280 FLOAT* WINAPI D3DXSHAdd(FLOAT *out, UINT order, const FLOAT *a, const FLOAT *b)
2281 {
2282     UINT i;
2283
2284     TRACE("out %p, order %u, a %p, b %p\n", out, order, a, b);
2285
2286     for (i = 0; i < order * order; i++)
2287         out[i] = a[i] + b[i];
2288
2289     return out;
2290 }
2291
2292 FLOAT WINAPI D3DXSHDot(UINT order, const FLOAT *a, const FLOAT *b)
2293 {
2294     FLOAT s;
2295     UINT i;
2296
2297     TRACE("order %u, a %p, b %p\n", order, a, b);
2298
2299     s = a[0] * b[0];
2300     for (i = 1; i < order * order; i++)
2301         s += a[i] * b[i];
2302
2303     return s;
2304 }
2305
2306 static void weightedcapintegrale(FLOAT *out, FLOAT order, FLOAT angle)
2307 {
2308     FLOAT coeff[3];
2309
2310     coeff[0] = cosf(angle);
2311
2312     out[0] = 2.0f * D3DX_PI * (1.0f - coeff[0]);
2313     out[1] = D3DX_PI * sinf(angle) * sinf(angle);
2314     if (order <= 2)
2315         return;
2316
2317     out[2] = coeff[0] * out[1];
2318     if (order == 3)
2319         return;
2320
2321     coeff[1] = coeff[0] * coeff[0];
2322     coeff[2] = coeff[1] * coeff[1];
2323
2324     out[3] = D3DX_PI * (-1.25f * coeff[2] + 1.5f * coeff[1] - 0.25f);
2325     if (order == 4)
2326         return;
2327
2328     out[4] = -0.25f * D3DX_PI * coeff[0] * (7.0f * coeff[2] - 10.0f * coeff[1] + 3.0f);
2329     if (order == 5)
2330         return;
2331
2332     out[5] = D3DX_PI * (-2.625f * coeff[2] * coeff[1] + 4.375f * coeff[2] - 1.875f * coeff[1] + 0.125f);
2333 }
2334
2335 HRESULT WINAPI D3DXSHEvalConeLight(UINT order, const D3DXVECTOR3 *dir, FLOAT radius,
2336     FLOAT Rintensity, FLOAT Gintensity, FLOAT Bintensity, FLOAT *rout, FLOAT *gout, FLOAT *bout)
2337 {
2338     FLOAT cap[6], clamped_angle, norm, scale, temp;
2339     UINT i, index, j;
2340
2341     TRACE("order %u, dir %p, radius %f, red %f, green %f, blue %f, rout %p, gout %p, bout %p\n",
2342         order, dir, radius, Rintensity, Gintensity, Bintensity, rout, gout, bout);
2343
2344     if (radius <= 0.0f)
2345         return D3DXSHEvalDirectionalLight(order, dir, Rintensity, Gintensity, Bintensity, rout, gout, bout);
2346
2347     clamped_angle = (radius > D3DX_PI / 2.0f) ? (D3DX_PI / 2.0f) : radius;
2348     norm = sinf(clamped_angle) * sinf(clamped_angle);
2349
2350     if (order > D3DXSH_MAXORDER)
2351     {
2352         WARN("Order clamped at D3DXSH_MAXORDER\n");
2353         order = D3DXSH_MAXORDER;
2354     }
2355
2356     weightedcapintegrale(cap, order, radius);
2357     D3DXSHEvalDirection(rout, order, dir);
2358
2359     for (i = 0; i < order; i++)
2360     {
2361         scale = cap[i] / norm;
2362
2363         for (j = 0; j < 2 * i + 1; j++)
2364         {
2365             index = i * i + j;
2366             temp = rout[index] * scale;
2367
2368             rout[index] = temp * Rintensity;
2369             if (gout)
2370                 gout[index] = temp * Gintensity;
2371             if (bout)
2372                 bout[index] = temp * Bintensity;
2373         }
2374     }
2375
2376     return D3D_OK;
2377 }
2378
2379 FLOAT* WINAPI D3DXSHEvalDirection(FLOAT *out, UINT order, const D3DXVECTOR3 *dir)
2380 {
2381     const FLOAT dirxx = dir->x * dir->x;
2382     const FLOAT dirxy = dir->x * dir->y;
2383     const FLOAT dirxz = dir->x * dir->z;
2384     const FLOAT diryy = dir->y * dir->y;
2385     const FLOAT diryz = dir->y * dir->z;
2386     const FLOAT dirzz = dir->z * dir->z;
2387     const FLOAT dirxxxx = dirxx * dirxx;
2388     const FLOAT diryyyy = diryy * diryy;
2389     const FLOAT dirzzzz = dirzz * dirzz;
2390     const FLOAT dirxyxy = dirxy * dirxy;
2391
2392     TRACE("out %p, order %u, dir %p\n", out, order, dir);
2393
2394     if ((order < D3DXSH_MINORDER) || (order > D3DXSH_MAXORDER))
2395         return out;
2396
2397     out[0] = 0.5f / sqrtf(D3DX_PI);
2398     out[1] = -0.5f / sqrtf(D3DX_PI / 3.0f) * dir->y;
2399     out[2] = 0.5f / sqrtf(D3DX_PI / 3.0f) * dir->z;
2400     out[3] = -0.5f / sqrtf(D3DX_PI / 3.0f) * dir->x;
2401     if (order == 2)
2402         return out;
2403
2404     out[4] = 0.5f / sqrtf(D3DX_PI / 15.0f) * dirxy;
2405     out[5] = -0.5f / sqrtf(D3DX_PI / 15.0f) * diryz;
2406     out[6] = 0.25f / sqrtf(D3DX_PI / 5.0f) * (3.0f * dirzz - 1.0f);
2407     out[7] = -0.5f / sqrtf(D3DX_PI / 15.0f) * dirxz;
2408     out[8] = 0.25f / sqrtf(D3DX_PI / 15.0f) * (dirxx - diryy);
2409     if (order == 3)
2410         return out;
2411
2412     out[9] = -sqrtf(70.0f / D3DX_PI) / 8.0f * dir->y * (3.0f * dirxx - diryy);
2413     out[10] = sqrtf(105.0f / D3DX_PI) / 2.0f * dirxy * dir->z;
2414     out[11] = -sqrtf(42.0 / D3DX_PI) / 8.0f * dir->y * (-1.0f + 5.0f * dirzz);
2415     out[12] = sqrtf(7.0f / D3DX_PI) / 4.0f * dir->z * (5.0f * dirzz - 3.0f);
2416     out[13] = sqrtf(42.0 / D3DX_PI) / 8.0f * dir->x * (1.0f - 5.0f * dirzz);
2417     out[14] = sqrtf(105.0f / D3DX_PI) / 4.0f * dir->z * (dirxx - diryy);
2418     out[15] = -sqrtf(70.0f / D3DX_PI) / 8.0f * dir->x * (dirxx - 3.0f * diryy);
2419     if (order == 4)
2420         return out;
2421
2422     out[16] = 0.75f * sqrtf(35.0f / D3DX_PI) * dirxy * (dirxx - diryy);
2423     out[17] = 3.0f * dir->z * out[9];
2424     out[18] = 0.75f * sqrtf(5.0f / D3DX_PI) * dirxy * (7.0f * dirzz - 1.0f);
2425     out[19] = 0.375f * sqrtf(10.0f / D3DX_PI) * diryz * (3.0f - 7.0f * dirzz);
2426     out[20] = 3.0f / (16.0f * sqrtf(D3DX_PI)) * (35.0f * dirzzzz - 30.f * dirzz + 3.0f);
2427     out[21] = 0.375f * sqrtf(10.0f / D3DX_PI) * dirxz * (3.0f - 7.0f * dirzz);
2428     out[22] = 0.375f * sqrtf(5.0f / D3DX_PI) * (dirxx - diryy) * (7.0f * dirzz - 1.0f);
2429     out[23] = 3.0 * dir->z * out[15];
2430     out[24] = 3.0f / 16.0f * sqrtf(35.0f / D3DX_PI) * (dirxxxx - 6.0f * dirxyxy + diryyyy);
2431     if (order == 5)
2432         return out;
2433
2434     out[25] = -3.0f/ 32.0f * sqrtf(154.0f / D3DX_PI) * dir->y * (5.0f * dirxxxx - 10.0f * dirxyxy + diryyyy);
2435     out[26] = 0.75f * sqrtf(385.0f / D3DX_PI) * dirxy * dir->z * (dirxx - diryy);
2436     out[27] = sqrtf(770.0f / D3DX_PI) / 32.0f * dir->y * (3.0f * dirxx - diryy) * (1.0f - 9.0f * dirzz);
2437     out[28] = sqrtf(1155.0f / D3DX_PI) / 4.0f * dirxy * dir->z * (3.0f * dirzz - 1.0f);
2438     out[29] = sqrtf(165.0f / D3DX_PI) / 16.0f * dir->y * (14.0f * dirzz - 21.0f * dirzzzz - 1.0f);
2439     out[30] = sqrtf(11.0f / D3DX_PI) / 16.0f * dir->z * (63.0f * dirzzzz - 70.0f * dirzz + 15.0f);
2440     out[31] = sqrtf(165.0f / D3DX_PI) / 16.0f * dir->x * (14.0f * dirzz - 21.0f * dirzzzz - 1.0f);
2441     out[32] = sqrtf(1155.0f / D3DX_PI) / 8.0f * dir->z * (dirxx - diryy) * (3.0f * dirzz - 1.0f);
2442     out[33] = sqrtf(770.0f / D3DX_PI) / 32.0f * dir->x * (dirxx - 3.0f * diryy) * (1.0f - 9.0f * dirzz);
2443     out[34] = 3.0f / 16.0f * sqrtf(385.0f / D3DX_PI) * dir->z * (dirxxxx - 6.0 * dirxyxy + diryyyy);
2444     out[35] = -3.0f/ 32.0f * sqrtf(154.0f / D3DX_PI) * dir->x * (dirxxxx - 10.0f * dirxyxy + 5.0f * diryyyy);
2445
2446     return out;
2447 }
2448
2449 HRESULT WINAPI D3DXSHEvalDirectionalLight(UINT order, const D3DXVECTOR3 *dir, FLOAT Rintensity, FLOAT Gintensity, FLOAT Bintensity, FLOAT *Rout, FLOAT *Gout, FLOAT *Bout)
2450 {
2451     FLOAT s, temp;
2452     UINT j;
2453
2454     TRACE("Order %u, Vector %p, Red %f, Green %f, Blue %f, Rout %p, Gout %p, Bout %p\n", order, dir, Rintensity, Gintensity, Bintensity, Rout, Gout, Bout);
2455
2456     s = 0.75f;
2457     if ( order > 2 )
2458         s += 5.0f / 16.0f;
2459     if ( order > 4 )
2460         s -= 3.0f / 32.0f;
2461     s /= D3DX_PI;
2462
2463     D3DXSHEvalDirection(Rout, order, dir);
2464     for (j = 0; j < order * order; j++)
2465     {
2466         temp = Rout[j] / s;
2467
2468         Rout[j] = Rintensity * temp;
2469         if ( Gout )
2470             Gout[j] = Gintensity * temp;
2471         if ( Bout )
2472             Bout[j] = Bintensity * temp;
2473     }
2474
2475     return D3D_OK;
2476 }
2477
2478 HRESULT WINAPI D3DXSHEvalHemisphereLight(UINT order, const D3DXVECTOR3 *dir, D3DXCOLOR top, D3DXCOLOR bottom,
2479     FLOAT *rout, FLOAT *gout, FLOAT *bout)
2480 {
2481     FLOAT a[2], temp[4];
2482     UINT i, j;
2483
2484     TRACE("order %u, dir %p, rout %p, gout %p, bout %p\n", order, dir, rout, gout, bout);
2485
2486     D3DXSHEvalDirection(temp, 2, dir);
2487
2488     a[0] = (top.r + bottom.r) * 3.0f * D3DX_PI;
2489     a[1] = (top.r - bottom.r) * D3DX_PI;
2490     for (i = 0; i < order; i++)
2491         for (j = 0; j < 2 * i + 1; j++)
2492             if (i < 2)
2493                 rout[i * i + j] = temp[i * i + j] * a[i];
2494             else
2495                 rout[i * i + j] = 0.0f;
2496
2497     if (gout)
2498     {
2499         a[0] = (top.g + bottom.g) * 3.0f * D3DX_PI;
2500         a[1] = (top.g - bottom.g) * D3DX_PI;
2501         for (i = 0; i < order; i++)
2502             for (j = 0; j < 2 * i + 1; j++)
2503                 if (i < 2)
2504                     gout[i * i + j] = temp[i * i + j] * a[i];
2505                 else
2506                     gout[i * i + j] = 0.0f;
2507     }
2508
2509     if (bout)
2510     {
2511         a[0] = (top.b + bottom.b) * 3.0f * D3DX_PI;
2512         a[1] = (top.b - bottom.b) * D3DX_PI;
2513         for (i = 0; i < order; i++)
2514             for (j = 0; j < 2 * i + 1; j++)
2515                 if (i < 2)
2516                     bout[i * i + j] = temp[i * i + j] * a[i];
2517                 else
2518                     bout[i * i + j] = 0.0f;
2519     }
2520
2521     return D3D_OK;
2522 }
2523
2524 HRESULT WINAPI D3DXSHEvalSphericalLight(UINT order, const D3DXVECTOR3 *dir, FLOAT radius,
2525     FLOAT Rintensity, FLOAT Gintensity, FLOAT Bintensity, FLOAT *rout, FLOAT *gout, FLOAT *bout)
2526 {
2527     D3DXVECTOR3 normal;
2528     FLOAT cap[6], clamped_angle, dist, temp;
2529     UINT i, index, j;
2530
2531     TRACE("order %u, dir %p, radius %f, red %f, green %f, blue %f, rout %p, gout %p, bout %p\n",
2532         order, dir, radius, Rintensity, Gintensity, Bintensity, rout, gout, bout);
2533
2534     if (order > D3DXSH_MAXORDER)
2535     {
2536         WARN("Order clamped at D3DXSH_MAXORDER\n");
2537         order = D3DXSH_MAXORDER;
2538     }
2539
2540     if (radius < 0.0f)
2541         radius = -radius;
2542
2543     dist = D3DXVec3Length(dir);
2544     clamped_angle = (dist <= radius) ? D3DX_PI / 2.0f : asinf(radius / dist);
2545
2546     weightedcapintegrale(cap, order, clamped_angle);
2547     D3DXVec3Normalize(&normal, dir);
2548     D3DXSHEvalDirection(rout, order, &normal);
2549
2550     for (i = 0; i < order; i++)
2551         for (j = 0; j < 2 * i + 1; j++)
2552         {
2553             index = i * i + j;
2554             temp = rout[index] * cap[i];
2555
2556             rout[index] = temp * Rintensity;
2557             if (gout)
2558                 gout[index] = temp * Gintensity;
2559             if (bout)
2560                 bout[index] = temp * Bintensity;
2561         }
2562
2563     return D3D_OK;
2564 }
2565
2566 FLOAT * WINAPI D3DXSHMultiply2(FLOAT *out, const FLOAT *a, const FLOAT *b)
2567 {
2568     FLOAT ta, tb;
2569
2570     TRACE("out %p, a %p, b %p\n", out, a, b);
2571
2572     ta = 0.28209479f * a[0];
2573     tb = 0.28209479f * b[0];
2574
2575     out[0] = 0.28209479f * D3DXSHDot(2, a, b);
2576     out[1] = ta * b[1] + tb * a[1];
2577     out[2] = ta * b[2] + tb * a[2];
2578     out[3] = ta * b[3] + tb * a[3];
2579
2580     return out;
2581 }
2582
2583 FLOAT * WINAPI D3DXSHMultiply3(FLOAT *out, const FLOAT *a, const FLOAT *b)
2584 {
2585     FLOAT t, ta, tb;
2586
2587     TRACE("out %p, a %p, b %p\n", out, a, b);
2588
2589     out[0] = 0.28209479f * a[0] * b[0];
2590
2591     ta = 0.28209479f * a[0] - 0.12615662f * a[6] - 0.21850968f * a[8];
2592     tb = 0.28209479f * b[0] - 0.12615662f * b[6] - 0.21850968f * b[8];
2593     out[1] = ta * b[1] + tb * a[1];
2594     t = a[1] * b[1];
2595     out[0] += 0.28209479f * t;
2596     out[6] = -0.12615662f * t;
2597     out[8] = -0.21850968f * t;
2598
2599     ta = 0.21850968f * a[5];
2600     tb = 0.21850968f * b[5];
2601     out[1] += ta * b[2] + tb * a[2];
2602     out[2] = ta * b[1] + tb * a[1];
2603     t = a[1] * b[2] +a[2] * b[1];
2604     out[5] = 0.21850968f * t;
2605
2606     ta = 0.21850968f * a[4];
2607     tb = 0.21850968f * b[4];
2608     out[1] += ta * b[3] + tb * a[3];
2609     out[3]  = ta * b[1] + tb * a[1];
2610     t = a[1] * b[3] + a[3] * b[1];
2611     out[4] = 0.21850968f * t;
2612
2613     ta = 0.28209480f * a[0] + 0.25231326f * a[6];
2614     tb = 0.28209480f * b[0] + 0.25231326f * b[6];
2615     out[2] += ta * b[2] + tb * a[2];
2616     t = a[2] * b[2];
2617     out[0] += 0.28209480f * t;
2618     out[6] += 0.25231326f * t;
2619
2620     ta = 0.21850969f * a[7];
2621     tb = 0.21850969f * b[7];
2622     out[2] += ta * b[3] + tb * a[3];
2623     out[3] += ta * b[2] + tb * a[2];
2624     t = a[2] * b[3] + a[3] * b[2];
2625     out[7] = 0.21850969f * t;
2626
2627     ta = 0.28209479f * a[0] - 0.12615663f * a[6] + 0.21850969f * a[8];
2628     tb = 0.28209479f * b[0] - 0.12615663f * b[6] + 0.21850969f * b[8];
2629     out[3] += ta * b[3] + tb * a[3];
2630     t = a[3] * b[3];
2631     out[0] += 0.28209479f * t;
2632     out[6] -= 0.12615663f * t;
2633     out[8] += 0.21850969f * t;
2634
2635     ta = 0.28209479f * a[0] - 0.18022375f * a[6];
2636     tb = 0.28209479f * b[0] - 0.18022375f * b[6];
2637     out[4] += ta * b[4] + tb * a[4];
2638     t = a[4] * b[4];
2639     out[0] += 0.28209479f * t;
2640     out[6] -= 0.18022375f * t;
2641
2642     ta = 0.15607835f * a[7];
2643     tb = 0.15607835f * b[7];
2644     out[4] += ta * b[5] + tb * a[5];
2645     out[5] += ta * b[4] + tb * a[4];
2646     t = a[4] * b[5] + a[5] * b[4];
2647     out[7] += 0.15607834f * t;
2648
2649     ta = 0.28209479f * a[0] + 0.09011186 * a[6] - 0.15607835f * a[8];
2650     tb = 0.28209479f * b[0] + 0.09011186 * b[6] - 0.15607835f * b[8];
2651     out[5] += ta * b[5] + tb * a[5];
2652     t = a[5] * b[5];
2653     out[0] += 0.28209479f * t;
2654     out[6] += 0.09011186f * t;
2655     out[8] -= 0.15607835f * t;
2656
2657     ta = 0.28209480f * a[0];
2658     tb = 0.28209480f * b[0];
2659     out[6] += ta * b[6] + tb * a[6];
2660     t = a[6] * b[6];
2661     out[0] += 0.28209480f * t;
2662     out[6] += 0.18022376f * t;
2663
2664     ta = 0.28209479f * a[0] + 0.09011186 * a[6] + 0.15607835f * a[8];
2665     tb = 0.28209479f * b[0] + 0.09011186 * b[6] + 0.15607835f * b[8];
2666     out[7] += ta * b[7] + tb * a[7];
2667     t = a[7] * b[7];
2668     out[0] += 0.28209479f * t;
2669     out[6] += 0.09011186f * t;
2670     out[8] += 0.15607835f * t;
2671
2672     ta = 0.28209479f * a[0] - 0.18022375f * a[6];
2673     tb = 0.28209479f * b[0] - 0.18022375f * b[6];
2674     out[8] += ta * b[8] + tb * a[8];
2675     t = a[8] * b[8];
2676     out[0] += 0.28209479f * t;
2677     out[6] -= 0.18022375f * t;
2678
2679     return out;
2680 }
2681
2682 FLOAT * WINAPI D3DXSHMultiply4(FLOAT *out, const FLOAT *a, const FLOAT *b)
2683 {
2684     FLOAT ta, tb, t;
2685
2686     TRACE("out %p, a %p, b %p\n", out, a, b);
2687
2688     out[0] = 0.28209479f * a[0] * b[0];
2689
2690     ta = 0.28209479f * a[0] - 0.12615663f * a[6] - 0.21850969f * a[8];
2691     tb = 0.28209479f * b[0] - 0.12615663f * b[6] - 0.21850969f * b[8];
2692     out[1] = ta * b[1] + tb * a[1];
2693     t = a[1] * b[1];
2694     out[0] += 0.28209479f * t;
2695     out[6] = -0.12615663f * t;
2696     out[8] = -0.21850969f * t;
2697
2698     ta = 0.21850969f * a[3] - 0.05839917f * a[13] - 0.22617901f * a[15];
2699     tb = 0.21850969f * b[3] - 0.05839917f * b[13] - 0.22617901f * b[15];
2700     out[1] += ta * b[4] + tb * a[4];
2701     out[4] = ta * b[1] + tb * a[1];
2702     t = a[1] * b[4] + a[4] * b[1];
2703     out[3] = 0.21850969f * t;
2704     out[13] = -0.05839917f * t;
2705     out[15] = -0.22617901f * t;
2706
2707     ta = 0.21850969f * a[2] - 0.14304817f * a[12] - 0.18467439f * a[14];
2708     tb = 0.21850969f * b[2] - 0.14304817f * b[12] - 0.18467439f * b[14];
2709     out[1] += ta * b[5] + tb * a[5];
2710     out[5] = ta * b[1] + tb * a[1];
2711     t = a[1] * b[5] + a[5] * b[1];
2712     out[2] = 0.21850969f * t;
2713     out[12] = -0.14304817f * t;
2714     out[14] = -0.18467439f * t;
2715
2716     ta = 0.20230066f * a[11];
2717     tb = 0.20230066f * b[11];
2718     out[1] += ta * b[6] + tb * a[6];
2719     out[6] += ta * b[1] + tb * a[1];
2720     t = a[1] * b[6] + a[6] * b[1];
2721     out[11] = 0.20230066f * t;
2722
2723     ta = 0.22617901f * a[9] + 0.05839917f * a[11];
2724     tb = 0.22617901f * b[9] + 0.05839917f * b[11];
2725     out[1] += ta * b[8] + tb * a[8];
2726     out[8] += ta * b[1] + tb * a[1];
2727     t = a[1] * b[8] + a[8] * b[1];
2728     out[9] = 0.22617901f * t;
2729     out[11] += 0.05839917f * t;
2730
2731     ta = 0.28209480f * a[0] + 0.25231326f * a[6];
2732     tb = 0.28209480f * b[0] + 0.25231326f * b[6];
2733     out[2] += ta * b[2] + tb * a[2];
2734     t = a[2] * b[2];
2735     out[0] += 0.28209480f * t;
2736     out[6] += 0.25231326f * t;
2737
2738     ta = 0.24776671f * a[12];
2739     tb = 0.24776671f * b[12];
2740     out[2] += ta * b[6] + tb * a[6];
2741     out[6] += ta * b[2] + tb * a[2];
2742     t = a[2] * b[6] + a[6] * b[2];
2743     out[12] += 0.24776671f * t;
2744
2745     ta = 0.28209480f * a[0] - 0.12615663f * a[6] + 0.21850969f * a[8];
2746     tb = 0.28209480f * b[0] - 0.12615663f * b[6] + 0.21850969f * b[8];
2747     out[3] += ta * b[3] + tb * a[3];
2748     t = a[3] * b[3];
2749     out[0] += 0.28209480f * t;
2750     out[6] -= 0.12615663f * t;
2751     out[8] += 0.21850969f * t;
2752
2753     ta = 0.20230066f * a[13];
2754     tb = 0.20230066f * b[13];
2755     out[3] += ta * b[6] + tb * a[6];
2756     out[6] += ta * b[3] + tb * a[3];
2757     t = a[3] * b[6] + a[6] * b[3];
2758     out[13] += 0.20230066f * t;
2759
2760     ta = 0.21850969f * a[2] - 0.14304817f * a[12] + 0.18467439f * a[14];
2761     tb = 0.21850969f * b[2] - 0.14304817f * b[12] + 0.18467439f * b[14];
2762     out[3] += ta * b[7] + tb * a[7];
2763     out[7] = ta * b[3] + tb * a[3];
2764     t = a[3] * b[7] + a[7] * b[3];
2765     out[2] += 0.21850969f * t;
2766     out[12] -= 0.14304817f * t;
2767     out[14] += 0.18467439f * t;
2768
2769     ta = -0.05839917f * a[13] + 0.22617901f * a[15];
2770     tb = -0.05839917f * b[13] + 0.22617901f * b[15];
2771     out[3] += ta * b[8] + tb * a[8];
2772     out[8] += ta * b[3] + tb * a[3];
2773     t = a[3] * b[8] + a[8] * b[3];
2774     out[13] -= 0.05839917f * t;
2775     out[15] += 0.22617901f * t;
2776
2777     ta = 0.28209479f * a[0] - 0.18022375f * a[6];
2778     tb = 0.28209479f * b[0] - 0.18022375f * b[6];
2779     out[4] += ta * b[4] + tb * a[4];
2780     t = a[4] * b[4];
2781     out[0] += 0.28209479f * t;
2782     out[6] -= 0.18022375f * t;
2783
2784     ta = 0.15607835f * a[7];
2785     tb = 0.15607835f * b[7];
2786     out[4] += ta * b[5] + tb * a[5];
2787     out[5] += ta * b[4] + tb * a[4];
2788     t = a[4] * b[5] + a[5] * b[4];
2789     out[7] += 0.15607835f * t;
2790
2791     ta = 0.22617901f * a[3] - 0.09403160f * a[13];
2792     tb = 0.22617901f * b[3] - 0.09403160f * b[13];
2793     out[4] += ta * b[9] + tb * a[9];
2794     out[9] += ta * b[4] + tb * a[4];
2795     t = a[4] * b[9] + a[9] * b[4];
2796     out[3] += 0.22617901f * t;
2797     out[13] -= 0.09403160f * t;
2798
2799     ta = 0.18467439f * a[2] - 0.18806319f * a[12];
2800     tb = 0.18467439f * b[2] - 0.18806319f * b[12];
2801     out[4] += ta * b[10] + tb * a [10];
2802     out[10] = ta * b[4] + tb * a[4];
2803     t = a[4] * b[10] + a[10] * b[4];
2804     out[2] += 0.18467439f * t;
2805     out[12] -= 0.18806319f * t;
2806
2807     ta = -0.05839917f * a[3] + 0.14567312f * a[13] + 0.09403160f * a[15];
2808     tb = -0.05839917f * b[3] + 0.14567312f * b[13] + 0.09403160f * b[15];
2809     out[4] += ta * b[11] + tb * a[11];
2810     out[11] += ta * b[4] + tb * a[4];
2811     t = a[4] * b[11] + a[11] * b[4];
2812     out[3] -= 0.05839917f * t;
2813     out[13] += 0.14567312f * t;
2814     out[15] += 0.09403160f * t;
2815
2816     ta = 0.28209479f * a[0] + 0.09011186f * a[6] - 0.15607835f * a[8];
2817     tb = 0.28209479f * b[0] + 0.09011186f * b[6] - 0.15607835f * b[8];
2818     out[5] += ta * b[5] + tb * a[5];
2819     t = a[5] * b[5];
2820     out[0] += 0.28209479f * t;
2821     out[6] += 0.09011186f * t;
2822     out[8] -= 0.15607835f * t;
2823
2824     ta = 0.14867701f * a[14];
2825     tb = 0.14867701f * b[14];
2826     out[5] += ta * b[9] + tb * a[9];
2827     out[9] += ta * b[5] + tb * a[5];
2828     t = a[5] * b[9] + a[9] * b[5];
2829     out[14] += 0.14867701f * t;
2830
2831     ta = 0.18467439f * a[3] + 0.11516472f * a[13] - 0.14867701f * a[15];
2832     tb = 0.18467439f * b[3] + 0.11516472f * b[13] - 0.14867701f * b[15];
2833     out[5] += ta * b[10] + tb * a[10];
2834     out[10] += ta * b[5] + tb * a[5];
2835     t = a[5] * b[10] + a[10] * b[5];
2836     out[3] += 0.18467439f * t;
2837     out[13] += 0.11516472f * t;
2838     out[15] -= 0.14867701f * t;
2839
2840     ta = 0.23359668f * a[2] + 0.05947080f * a[12] - 0.11516472f * a[14];
2841     tb = 0.23359668f * b[2] + 0.05947080f * b[12] - 0.11516472f * b[14];
2842     out[5] += ta * b[11] + tb * a[11];
2843     out[11] += ta * b[5] + tb * a[5];
2844     t = a[5] * b[11] + a[11] * b[5];
2845     out[2] += 0.23359668f * t;
2846     out[12] += 0.05947080f * t;
2847     out[14] -= 0.11516472f * t;
2848
2849     ta = 0.28209479f * a[0];
2850     tb = 0.28209479f * b[0];
2851     out[6] += ta * b[6] + tb * a[6];
2852     t = a[6] * b[6];
2853     out[0] += 0.28209479f * t;
2854     out[6] += 0.18022376f * t;
2855
2856     ta = 0.09011186f * a[6] + 0.28209479f * a[0] + 0.15607835f * a[8];
2857     tb = 0.09011186f * b[6] + 0.28209479f * b[0] + 0.15607835f * b[8];
2858     out[7] += ta * b[7] + tb * a[7];
2859     t = a[7] * b[7];
2860     out[6] += 0.09011186f * t;
2861     out[0] += 0.28209479f * t;
2862     out[8] += 0.15607835f * t;
2863
2864     ta = 0.14867701f * a[9] + 0.18467439f * a[1] + 0.11516472f * a[11];
2865     tb = 0.14867701f * b[9] + 0.18467439f * b[1] + 0.11516472f * b[11];
2866     out[7] += ta * b[10] + tb * a[10];
2867     out[10] += ta * b[7] + tb * a[7];
2868     t = a[7] * b[10] + a[10] * b[7];
2869     out[9] += 0.14867701f * t;
2870     out[1] += 0.18467439f * t;
2871     out[11] += 0.11516472f * t;
2872
2873     ta = 0.05947080f * a[12] + 0.23359668f * a[2] + 0.11516472f * a[14];
2874     tb = 0.05947080f * b[12] + 0.23359668f * b[2] + 0.11516472f * b[14];
2875     out[7] += ta * b[13] + tb * a[13];
2876     out[13] += ta * b[7]+ tb * a[7];
2877     t = a[7] * b[13] + a[13] * b[7];
2878     out[12] += 0.05947080f * t;
2879     out[2] += 0.23359668f * t;
2880     out[14] += 0.11516472f * t;
2881
2882     ta = 0.14867701f * a[15];
2883     tb = 0.14867701f * b[15];
2884     out[7] += ta * b[14] + tb * a[14];
2885     out[14] += ta * b[7] + tb * a[7];
2886     t = a[7] * b[14] + a[14] * b[7];
2887     out[15] += 0.14867701f * t;
2888
2889     ta = 0.28209479f * a[0] - 0.18022375f * a[6];
2890     tb = 0.28209479f * b[0] - 0.18022375f * b[6];
2891     out[8] += ta * b[8] + tb * a[8];
2892     t = a[8] * b[8];
2893     out[0] += 0.28209479f * t;
2894     out[6] -= 0.18022375f * t;
2895
2896     ta = -0.09403160f * a[11];
2897     tb = -0.09403160f * b[11];
2898     out[8] += ta * b[9] + tb * a[9];
2899     out[9] += ta * b[8] + tb * a[8];
2900     t = a[8] * b[9] + a[9] * b[8];
2901     out[11] -= 0.09403160f * t;
2902
2903     ta = -0.09403160f * a[15];
2904     tb = -0.09403160f * b[15];
2905     out[8] += ta * b[13] + tb * a[13];
2906     out[13] += ta * b[8] + tb * a[8];
2907     t = a[8] * b[13] + a[13] * b[8];
2908     out[15] -= 0.09403160f * t;
2909
2910     ta = 0.18467439f * a[2] - 0.18806319f * a[12];
2911     tb = 0.18467439f * b[2] - 0.18806319f * b[12];
2912     out[8] += ta * b[14] + tb * a[14];
2913     out[14] += ta * b[8] + tb * a[8];
2914     t = a[8] * b[14] + a[14] * b[8];
2915     out[2] += 0.18467439f * t;
2916     out[12] -= 0.18806319f * t;
2917
2918     ta = -0.21026104f * a[6] + 0.28209479f * a[0];
2919     tb = -0.21026104f * b[6] + 0.28209479f * b[0];
2920     out[9] += ta * b[9] + tb * a[9];
2921     t = a[9] * b[9];
2922     out[6] -= 0.21026104f * t;
2923     out[0] += 0.28209479f * t;
2924
2925     ta = 0.28209479f * a[0];
2926     tb = 0.28209479f * b[0];
2927     out[10] += ta * b[10] + tb * a[10];
2928     t = a[10] * b[10];
2929     out[0] += 0.28209479f * t;
2930
2931     ta = 0.28209479f * a[0] + 0.12615663f * a[6] - 0.14567312f * a[8];
2932     tb = 0.28209479f * b[0] + 0.12615663f * b[6] - 0.14567312f * b[8];
2933     out[11] += ta * b[11] + tb * a[11];
2934     t = a[11] * b[11];
2935     out[0] += 0.28209479f * t;
2936     out[6] += 0.12615663f * t;
2937     out[8] -= 0.14567312f * t;
2938
2939     ta = 0.28209479f * a[0] + 0.16820885f * a[6];
2940     tb = 0.28209479f * b[0] + 0.16820885f * b[6];
2941     out[12] += ta * b[12] + tb * a[12];
2942     t = a[12] * b[12];
2943     out[0] += 0.28209479f * t;
2944     out[6] += 0.16820885f * t;
2945
2946     ta =0.28209479f * a[0] + 0.14567312f * a[8] + 0.12615663f * a[6];
2947     tb =0.28209479f * b[0] + 0.14567312f * b[8] + 0.12615663f * b[6];
2948     out[13] += ta * b[13] + tb * a[13];
2949     t = a[13] * b[13];
2950     out[0] += 0.28209479f * t;
2951     out[8] += 0.14567312f * t;
2952     out[6] += 0.12615663f * t;
2953
2954     ta = 0.28209479f * a[0];
2955     tb = 0.28209479f * b[0];
2956     out[14] += ta * b[14] + tb * a[14];
2957     t = a[14] * b[14];
2958     out[0] += 0.28209479f * t;
2959
2960     ta = 0.28209479f * a[0] - 0.21026104f * a[6];
2961     tb = 0.28209479f * b[0] - 0.21026104f * b[6];
2962     out[15] += ta * b[15] + tb * a[15];
2963     t = a[15] * b[15];
2964     out[0] += 0.28209479f * t;
2965     out[6] -= 0.21026104f * t;
2966
2967     return out;
2968 }
2969
2970 static void rotate_X(FLOAT *out, UINT order, FLOAT a, FLOAT *in)
2971 {
2972     out[0] = in[0];
2973
2974     out[1] = a * in[2];
2975     out[2] = -a * in[1];
2976     out[3] = in[3];
2977
2978     out[4] = a * in[7];
2979     out[5] = -in[5];
2980     out[6] = -0.5f * in[6] - 0.8660253882f * in[8];
2981     out[7] = -a * in[4];
2982     out[8] = -0.8660253882f * in[6] + 0.5f * in[8];
2983     out[9] = -a * 0.7905694842f * in[12] + a * 0.6123724580f * in[14];
2984
2985     out[10] = -in[10];
2986     out[11] = -a * 0.6123724580f * in[12] - a * 0.7905694842f * in[14];
2987     out[12] = a * 0.7905694842f * in[9] + a * 0.6123724580f * in[11];
2988     out[13] = -0.25f * in[13] - 0.9682458639f * in[15];
2989     out[14] = -a * 0.6123724580f * in[9] + a * 0.7905694842f * in[11];
2990     out[15] = -0.9682458639f * in[13] + 0.25f * in[15];
2991     if (order == 4)
2992         return;
2993
2994     out[16] = -a * 0.9354143739f * in[21] + a * 0.3535533845f * in[23];
2995     out[17] = -0.75f * in[17] + 0.6614378095f * in[19];
2996     out[18] = -a * 0.3535533845f * in[21] - a * 0.9354143739f * in[23];
2997     out[19] = 0.6614378095f * in[17] + 0.75f * in[19];
2998     out[20] = 0.375f * in[20] + 0.5590170026f * in[22] + 0.7395099998f * in[24];
2999     out[21] = a * 0.9354143739f * in[16] + a * 0.3535533845f * in[18];
3000     out[22] = 0.5590170026f * in[20] + 0.5f * in[22] - 0.6614378691f * in[24];
3001     out[23] = -a * 0.3535533845f * in[16] + a * 0.9354143739f * in[18];
3002     out[24] = 0.7395099998f * in[20] - 0.6614378691f * in[22] + 0.125f * in[24];
3003     if (order == 5)
3004         return;
3005
3006     out[25] = a * 0.7015607357f * in[30] - a * 0.6846531630f * in[32] + a * 0.1976423711f * in[34];
3007     out[26] = -0.5f * in[26] + 0.8660253882f * in[28];
3008     out[27] = a * 0.5229125023f * in[30] + a * 0.3061861992f * in[32] - a * 0.7954951525 * in[34];
3009     out[28] = 0.8660253882f * in[26] + 0.5f * in[28];
3010     out[29] = a * 0.4841229022f * in[30] + a * 0.6614378691f * in[32] + a * 0.5728219748f * in[34];
3011     out[30] = -a * 0.7015607357f * in[25] - a * 0.5229125023f * in[27] - a * 0.4841229022f * in[29];
3012     out[31] = 0.125f * in[31] + 0.4050463140f * in[33] + 0.9057110548f * in[35];
3013     out[32] = a * 0.6846531630f * in[25] - a * 0.3061861992f * in[27] - a * 0.6614378691f * in[29];
3014     out[33] = 0.4050463140f * in[31] + 0.8125f * in[33] - 0.4192627370f * in[35];
3015     out[34] = -a * 0.1976423711f * in[25] + a * 0.7954951525f * in[27] - a * 0.5728219748f * in[29];
3016     out[35] = 0.9057110548f * in[31] - 0.4192627370f * in[33] + 0.0624999329f * in[35];
3017 }
3018
3019 FLOAT* WINAPI D3DXSHRotate(FLOAT *out, UINT order, const D3DXMATRIX *matrix, const FLOAT *in)
3020 {
3021     FLOAT alpha, beta, gamma, sinb, temp[36], temp1[36];
3022
3023     TRACE("out %p, order %u, matrix %p, in %p\n", out, order, matrix, in);
3024
3025     out[0] = in[0];
3026
3027     if ((order > D3DXSH_MAXORDER) || (order < D3DXSH_MINORDER))
3028         return out;
3029
3030     if (order <= 3)
3031     {
3032         out[1] = matrix->u.m[1][1] * in[1] - matrix->u.m[2][1] * in[2] + matrix->u.m[0][1] * in[3];
3033         out[2] = -matrix->u.m[1][2] * in[1] + matrix->u.m[2][2] * in[2] - matrix->u.m[0][2] * in[3];
3034         out[3] = matrix->u.m[1][0] * in[1] - matrix->u.m[2][0] * in[2] + matrix->u.m[0][0] * in[3];
3035
3036         if (order == 3)
3037         {
3038             FLOAT coeff[]={
3039                 matrix->u.m[1][0] * matrix->u.m[0][0], matrix->u.m[1][1] * matrix->u.m[0][1],
3040                 matrix->u.m[1][1] * matrix->u.m[2][1], matrix->u.m[1][0] * matrix->u.m[2][0],
3041                 matrix->u.m[2][0] * matrix->u.m[2][0], matrix->u.m[2][1] * matrix->u.m[2][1],
3042                 matrix->u.m[0][0] * matrix->u.m[2][0], matrix->u.m[0][1] * matrix->u.m[2][1],
3043                 matrix->u.m[0][1] * matrix->u.m[0][1], matrix->u.m[1][0] * matrix->u.m[1][0],
3044                 matrix->u.m[1][1] * matrix->u.m[1][1], matrix->u.m[0][0] * matrix->u.m[0][0], };
3045
3046             out[4] = (matrix->u.m[1][1] * matrix->u.m[0][0] + matrix->u.m[0][1] * matrix->u.m[1][0]) * in[4];
3047             out[4] -= (matrix->u.m[1][0] * matrix->u.m[2][1] + matrix->u.m[1][1] * matrix->u.m[2][0]) * in[5];
3048             out[4] += 1.7320508076f * matrix->u.m[2][0] * matrix->u.m[2][1] * in[6];
3049             out[4] -= (matrix->u.m[0][1] * matrix->u.m[2][0] + matrix->u.m[0][0] * matrix->u.m[2][1]) * in[7];
3050             out[4] += (matrix->u.m[0][0] * matrix->u.m[0][1] - matrix->u.m[1][0] * matrix->u.m[1][1]) * in[8];
3051
3052             out[5] = (matrix->u.m[1][1] * matrix->u.m[2][2] + matrix->u.m[1][2] * matrix->u.m[2][1]) * in[5];
3053             out[5] -= (matrix->u.m[1][1] * matrix->u.m[0][2] + matrix->u.m[1][2] * matrix->u.m[0][1]) * in[4];
3054             out[5] -= 1.7320508076f * matrix->u.m[2][2] * matrix->u.m[2][1] * in[6];
3055             out[5] += (matrix->u.m[0][2] * matrix->u.m[2][1] + matrix->u.m[0][1] * matrix->u.m[2][2]) * in[7];
3056             out[5] -= (matrix->u.m[0][1] * matrix->u.m[0][2] - matrix->u.m[1][1] * matrix->u.m[1][2]) * in[8];
3057
3058             out[6] = (matrix->u.m[2][2] * matrix->u.m[2][2] - 0.5f * (coeff[4] + coeff[5])) * in[6];
3059             out[6] -= (0.5773502692f * (coeff[0] + coeff[1]) - 1.1547005384f * matrix->u.m[1][2] * matrix->u.m[0][2]) * in[4];
3060             out[6] += (0.5773502692f * (coeff[2] + coeff[3]) - 1.1547005384f * matrix->u.m[1][2] * matrix->u.m[2][2]) * in[5];
3061             out[6] += (0.5773502692f * (coeff[6] + coeff[7]) - 1.1547005384f * matrix->u.m[0][2] * matrix->u.m[2][2]) * in[7];
3062             out[6] += (0.2886751347f * (coeff[9] - coeff[8] + coeff[10] - coeff[11]) - 0.5773502692f *
3063                   (matrix->u.m[1][2] * matrix->u.m[1][2] - matrix->u.m[0][2] * matrix->u.m[0][2])) * in[8];
3064
3065             out[7] = (matrix->u.m[0][0] * matrix->u.m[2][2] + matrix->u.m[0][2] * matrix->u.m[2][0]) * in[7];
3066             out[7] -= (matrix->u.m[1][0] * matrix->u.m[0][2] + matrix->u.m[1][2] * matrix->u.m[0][0]) * in[4];
3067             out[7] += (matrix->u.m[1][0] * matrix->u.m[2][2] + matrix->u.m[1][2] * matrix->u.m[2][0]) * in[5];
3068             out[7] -= 1.7320508076f * matrix->u.m[2][2] * matrix->u.m[2][0] * in[6];
3069             out[7] -= (matrix->u.m[0][0] * matrix->u.m[0][2] - matrix->u.m[1][0] * matrix->u.m[1][2]) * in[8];
3070
3071             out[8] = 0.5f * (coeff[11] - coeff[8] - coeff[9] + coeff[10]) * in[8];
3072             out[8] += (coeff[0] - coeff[1]) * in[4];
3073             out[8] += (coeff[2] - coeff[3]) * in[5];
3074             out[8] += 0.86602540f * (coeff[4] - coeff[5]) * in[6];
3075             out[8] += (coeff[7] - coeff[6]) * in[7];
3076         }
3077
3078         return out;
3079     }
3080
3081     if (fabsf(matrix->u.m[2][2]) != 1.0f)
3082     {
3083         sinb = sqrtf(1.0f - matrix->u.m[2][2] * matrix->u.m[2][2]);
3084         alpha = atan2f(matrix->u.m[2][1] / sinb, matrix->u.m[2][0] / sinb);
3085         beta = atan2f(sinb, matrix->u.m[2][2]);
3086         gamma = atan2f(matrix->u.m[1][2] / sinb, -matrix->u.m[0][2] / sinb);
3087     }
3088     else
3089     {
3090         alpha = atan2f(matrix->u.m[0][1], matrix->u.m[0][0]);
3091         beta = 0.0f;
3092         gamma = 0.0f;
3093     }
3094
3095     D3DXSHRotateZ(temp, order, gamma, in);
3096     rotate_X(temp1, order, 1.0f, temp);
3097     D3DXSHRotateZ(temp, order, beta, temp1);
3098     rotate_X(temp1, order, -1.0f, temp);
3099     D3DXSHRotateZ(out, order, alpha, temp1);
3100
3101     return out;
3102 }
3103
3104 FLOAT * WINAPI D3DXSHRotateZ(FLOAT *out, UINT order, FLOAT angle, const FLOAT *in)
3105 {
3106     UINT i, sum = 0;
3107     FLOAT c[5], s[5];
3108
3109     TRACE("out %p, order %u, angle %f, in %p\n", out, order, angle, in);
3110
3111     order = min(max(order, D3DXSH_MINORDER), D3DXSH_MAXORDER);
3112
3113     out[0] = in[0];
3114
3115     for (i = 1; i < order; i++)
3116     {
3117         UINT j;
3118
3119         c[i - 1] = cosf(i * angle);
3120         s[i - 1] = sinf(i * angle);
3121         sum += i * 2;
3122
3123         out[sum - i] = c[i - 1] * in[sum - i];
3124         out[sum - i] += s[i - 1] * in[sum + i];
3125         for (j = i - 1; j > 0; j--)
3126         {
3127             out[sum - j] = 0.0f;
3128             out[sum - j] = c[j - 1] * in[sum - j];
3129             out[sum - j] += s[j - 1] * in[sum + j];
3130         }
3131
3132         if (in == out)
3133             out[sum] = 0.0f;
3134         else
3135             out[sum] = in[sum];
3136
3137         for (j = 1; j < i; j++)
3138         {
3139             out[sum + j] = 0.0f;
3140             out[sum + j] = -s[j - 1] * in[sum - j];
3141             out[sum + j] += c[j - 1] * in[sum + j];
3142         }
3143         out[sum + i] = -s[i - 1] * in[sum - i];
3144         out[sum + i] += c[i - 1] * in[sum + i];
3145     }
3146
3147     return out;
3148 }
3149
3150 FLOAT* WINAPI D3DXSHScale(FLOAT *out, UINT order, const FLOAT *a, const FLOAT scale)
3151 {
3152     UINT i;
3153
3154     TRACE("out %p, order %u, a %p, scale %f\n", out, order, a, scale);
3155
3156     for (i = 0; i < order * order; i++)
3157         out[i] = a[i] * scale;
3158
3159     return out;
3160 }