Uniform discritization of SO(3) via random rotation matrices #28
Conversation
…into itr_ref_grid_so3
| rots = np.ones((n_rotations, 3, 3)) | ||
| rots = np.zeros((n_rotations, 3, 3)) | ||
| for i in range(n_rotations): | ||
| phi, theta, psi = np.random.random(3) * 2 * np.pi |
There was a problem hiding this comment.
This does not generate uniform rotations. They might look uniform in the picture, but they bunch at the poll. If you plotted a histogram on the sphere (like sky scrapers on the face of the globe) they would be higher on the polls.
Uniform is defined not in the parameters of the encoding of SO(3). I like to think of uniform on SO(3) as being defined as the distance between rotations tending to be equal. This raises the question of what distance is in SO(3). See the equation in "1. Introduction: Single Rotation Averaging" in Hartley, R., Trumpf, J., Dai, Y., & Li, H. (2013). Rotation Averaging. International Journal of Computer Vision, 103(3), 267–305. http://doi.org/10.1007/s11263-012-0601-0
You can use random_uniform in geomstats
https://geomstats.github.io/_modules/geomstats/geometry/special_orthogonal.html
I did something here that is uniform, but you should be able to just use what's in geomstats, not implement anything new yourselves. The less you implement, the less you have to test!
https://github.com/geoffwoollard/compSPI/blob/stash_simulate/src/coords.py#L61
There was a problem hiding this comment.
The last projection against (1,0,0) looks bunches at polls.
There was a problem hiding this comment.
Fixed distributions using geomstats
|
Oh good point I see the lint and test files are there in the github workflow folder... There just not in files changed since you are not merging into master, but back into |
|
Hmmm why isn't the linting/testing workflow being triggered. Will close and reopen. |
@geoffwoollard |
|
@thisTyler please fix the linting and testing issues. |
Codecov Report
@@ Coverage Diff @@
## dev #28 +/- ##
=========================================
Coverage 100.00% 100.00%
=========================================
Files 1 1
Lines 145 173 +28
=========================================
+ Hits 145 173 +28
Continue to review full report at Codecov.
|
| rots = np.ones((n_rotations, 3, 3)) | ||
| geom = special_orthogonal.SpecialOrthogonal(3, "matrix") | ||
| rots = geom.random_uniform(n_rotations) | ||
| negatives = np.tile(np.random.randint(2, size=n_rotations) * 2 - 1, (3, 3, 1)).T |
There was a problem hiding this comment.
I think you are negating every second rotation. A bit opaque in how the tiling works, so I would put a note in the docstring, or the variable name indicating how you are scaling each rotation by +1 or -1, with one random number per rotation.
There was a problem hiding this comment.
I'll add a note to the docstring. The reason for this is that geomstats only rotates on a half-sphere. So, I invert the rotations randomly in order to cover the entire sphere (I validated the functionality in a jupyter notebook and am happy to post images confirming this behaviour)
There was a problem hiding this comment.
Yes I recall this. Great work! Put as much as you can of that in the docstring and tests so that it doesn't go to waste.
The test has to be better than checking for the shape. For sure orthonormality
- orthogonality: each R.dot(R.T) is close np.eye(3)
- normality: np.linalg.det(R) is close 1
And I think something for uniformity would be great. If you can come up with some heuristic, would be a nice start. Doesn't have to be rigorous, just protest against a bad fail case. Some ideas (that work with large # of rotations)
- project to the three coordinates used to plot on sphere. Check roughly same amount above/below, through the three axes (up/down; left/right; front/back)
- convert to some encoding (Euler angles???) and check stats on those encoding (I think two Euler angles (ZYZ or ZXZ) are uniform and one is not)
- integrate with an even function (e.g. x2 + y2 + z2) that should be double the integration on half. If the function is odd (e.g. x3 + y3 + z3) it should integrate to roughly one.
- even: (axis_3_vectors[half_idx]**2).sum() is close 2*axis_3_vectors.sum()
- odd: (axis_3_vectors[half_idx]**3).sum() is close 0
Made issue here for ref: #39
|
I'll merge since everything is passing tests. Better tests can live on in issue #39 |
See #13.
Gridding of SO(3) via randomly distributed rotations - a 3D rotation matrix is generated from geomstats..
Demonstration of functionality (in each image 1000 rotation matrices are applied to a cartesian axis unit vector and plotted).
