lkilcher / dolfyn

A library for oceanographic doppler instruments such as Acoustic Doppler Profilers (ADPs, ADCPs) and Acoustic Doppler Velocimeters (ADVs).
BSD 3-Clause "New" or "Revised" License
42 stars 25 forks source link

User-input 'body2head_rotmat' read issue #65

Closed jmcvey3 closed 4 years ago

jmcvey3 commented 4 years ago

I was having an issue where changing the user-inputted 'body2head_rotmat' was making no difference within the motion correction, so I tracked it down to the private method _rotate_vel2body() (located in ./dolfyn/rotate/vector.py):

def _rotate_vel2body(advo):
    if (np.diag(np.eye(3)) == 1).all():
        advo.props['vel_rotated2body'] = True
    if 'vel_rotated2body' in advo.props and \
       advo.props['vel_rotated2body'] is True:
        # Don't re-rotate the data if its already been rotated.
        return
    if not rotb._check_rotmat_det(advo.props['body2head_rotmat']).all():
        raise ValueError("Invalid body-to-head rotation matrix"
                         " (determinant != 1).")
    # The transpose should do head to body.
    advo['vel'] = np.dot(advo.props['body2head_rotmat'].T, advo['vel'])
    advo.props['vel_rotated2body'] = True

The first "if" statement is never false, and subsequently the second "if" statement is always true, and the function returns without reaching the second to last line that does the head to body transpose.

I believe, for the first "if" statement, np.eye(3) was supposed to be advo.props['body2head_rotmat'], because if the latter = an identity matrix, the function appears to work fully as written.

So, replacing np.eye(3) with advo.props['body2head_rotmat'], along with using the right body2head orientation matrix appears to fix all the data that is run through dolfyn.adv.api.motion.correct_motion(). (It doesn't fix the uncorrected data - this private method to do the head-body transpose appears only to be called in dolfyn.adv.api.motion.correct_motion()?).

However this change does error out six of the nosetests because they can't find 'body2head_rotmat'; it doesn't appear to be loaded with the .h5 files? Not sure.

lkilcher commented 4 years ago

Wow! Really great job tracking this down! Thank you! I'll definitely try to fix this ASAP!

lkilcher commented 4 years ago

@jmcvey3 thanks again for flagging this! Note that I created #66 to solve this. Two questions tho:

lkilcher commented 4 years ago

Regarding my first question: can you send me the error trace you're getting?

jmcvey3 commented 4 years ago

Thanks, glad it helps! Regarding the second bullet, I wanted to overlay the corrected, uncorrected, and imu motion spectra all on the same plot. The coordinate rotations from instrument to earth to principal frames of reference don't work properly without that inital head-to-body rotation, since dolfyn defines the instrument frame as the ADV body frame.

## Perform motion correction
avm.motion.correct_motion(data, accel_filter, to_earth=True) # also rotates to earth frame

# Rotate the uncorrected data into the earth frame for comparison to motion correction:
data_despiked = data_despiked.rotate2('earth') # does not do head-to-body transpose

# Calc principal heading (from earth coordinates)
data.props['principal_heading'] = dlfn.calc_principal_heading(data.vel)
# use the corrected data's principal heading since it's the correct one
data_despiked.props['principal_heading'] = data.props['principal_heading']

# Then rotate it into a 'principal axes frame':
data = data.rotate2('principal')
data_despiked = data_despiked.rotate2('principal')

So I added _rotate_vel2body(advo) as the first step in inst2earth() method in the same Vector file (dolfyn/rotate/vector.py), so that when I rotate2('earth'), it does the head to body transpose automatically.

jmcvey3 commented 4 years ago

Not sure about the first bullet: I reinstalled and that fixed the errors. I must've changed something and forgotten about it. Though doing what I did in the above comment does error the nosetests because I don't think they're written to incorporate it.

lkilcher commented 4 years ago

@jmcvey3 I ended up making a much bigger change under #66 than my initial fix. Most of the extra work was in the category of putting stricter controls on what can/cannot be changed by the user in dat.props. This is something I'd been meaning to do for a while.

lkilcher commented 4 years ago

Resolved by #66.