I generally understand what it's intended to do ("Computes the inclination matrix I as well as the rotation matrix R transforming a vector from the device coordinate system to the world's coordinate system") but I don't get it how it works.
The javadoc is well written and the source is available here but I don't understand the math (for example, what is the mathematica/physical meaning of the values Hx, Hy and Hz? For example: Hx = Ey*Az - Ez*Ay). Also, what happens later in the method.
I've pasted the code from the GrepCode link above, leaving the source line numbers for easy reference. Thank you.
971 public static boolean More ...getRotationMatrix(float[] R, float[] I,
972 float[] gravity, float[] geomagnetic) {
973 // TODO: move this to native code for efficiency
974 float Ax = gravity[0];
975 float Ay = gravity[1];
976 float Az = gravity[2];
977 final float Ex = geomagnetic[0];
978 final float Ey = geomagnetic[1];
979 final float Ez = geomagnetic[2];
980 float Hx = Ey*Az - Ez*Ay;
981 float Hy = Ez*Ax - Ex*Az;
982 float Hz = Ex*Ay - Ey*Ax;
983 final float normH = (float)Math.sqrt(Hx*Hx + Hy*Hy + Hz*Hz);
984 if (normH < 0.1f) {
985 // device is close to free fall (or in space?), or close to
986 // magnetic north pole. Typical values are > 100.
987 return false;
988 }
989 final float invH = 1.0f / normH;
990 Hx *= invH;
991 Hy *= invH;
992 Hz *= invH;
993 final float invA = 1.0f / (float)Math.sqrt(Ax*Ax + Ay*Ay + Az*Az);
994 Ax *= invA;
995 Ay *= invA;
996 Az *= invA;
997 final float Mx = Ay*Hz - Az*Hy;
998 final float My = Az*Hx - Ax*Hz;
999 final float Mz = Ax*Hy - Ay*Hx;
1000 if (R != null) {
1001 if (R.length == 9) {
1002 R[0] = Hx; R[1] = Hy; R[2] = Hz;
1003 R[3] = Mx; R[4] = My; R[5] = Mz;
1004 R[6] = Ax; R[7] = Ay; R[8] = Az;
1005 } else if (R.length == 16) {
1006 R[0] = Hx; R[1] = Hy; R[2] = Hz; R[3] = 0;
1007 R[4] = Mx; R[5] = My; R[6] = Mz; R[7] = 0;
1008 R[8] = Ax; R[9] = Ay; R[10] = Az; R[11] = 0;
1009 R[12] = 0; R[13] = 0; R[14] = 0; R[15] = 1;
1010 }
1011 }
1012 if (I != null) {
1013 // compute the inclination matrix by projecting the geomagnetic
1014 // vector onto the Z (gravity) and X (horizontal component
1015 // of geomagnetic vector) axes.
1016 final float invE = 1.0f / (float)Math.sqrt(Ex*Ex + Ey*Ey + Ez*Ez);
1017 final float c = (Ex*Mx + Ey*My + Ez*Mz) * invE;
1018 final float s = (Ex*Ax + Ey*Ay + Ez*Az) * invE;
1019 if (I.length == 9) {
1020 I[0] = 1; I[1] = 0; I[2] = 0;
1021 I[3] = 0; I[4] = c; I[5] = s;
1022 I[6] = 0; I[7] =-s; I[8] = c;
1023 } else if (I.length == 16) {
1024 I[0] = 1; I[1] = 0; I[2] = 0;
1025 I[4] = 0; I[5] = c; I[6] = s;
1026 I[8] = 0; I[9] =-s; I[10]= c;
1027 I[3] = I[7] = I[11] = I[12] = I[13] = I[14] = 0;
1028 I[15] = 1;
1029 }
1030 }
1031 return true;
1032 }
They are the cross product. The assumption is one of the param passed in the method is gravity thus is the Sky axis. The other assumption is that the magnetic param passed in lies in the North-Sky plane. Thus the cross product of these two vectors is a vector orthogonal to the North-Sky plane which is East. Now the cross product of Sky and East is a vector orthogonal to both of these which is North. Normalize all of these will give a orthonormal basis for the World coordinate.
In the calculation above, H is East and M is North.