Moment refactor

This consolidates a lot of our previous attempts to compute moments. It implements all 3 iteration strategies, and includes a helper method to auto-choose a good strategy based on data size (this is currently very simplistic). I found and fixed a bunch of small bugs. 

There are a few open questions still:
  * Should moment 1 be returned in absolute or relative coordinates (ie separation from the first channel)? I've currently implemented things using absolute coordinates for velocity, and relative coordinates for RA/Dec. I think this is more formally correct per our discussion yesterday, but perhaps unfriendly for the case where the WCS is (nearly) cartesian and unrotated

* I'm using the equation `M(N) = Integral(Value * (pos - mom0) ** N dx) / Integral(Value dx)` for N >= 1. I'm not sure if it's more conveneint to return that to the power of 1/N (e.g. sigma instead of variance when N=2). Also, I don't know if that is the proper equation for N>2.

@keflavich @astrofrog @low-sky