Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

sigproc writer block needs some love #110

Open
telegraphic opened this issue Nov 15, 2017 · 3 comments
Open

sigproc writer block needs some love #110

telegraphic opened this issue Nov 15, 2017 · 3 comments

Comments

@telegraphic
Copy link
Collaborator

The current sigproc reader block a) doesn't treat angles correctly, and b) will crash in most cases! My fix so far has been to use astropy, but that adds a dependency, so starting an issue to discuss merits.

Angles

The sigproc definition for angles is actually a float, HHMMSS.ss. For example, 17h45m00.12s needs to be written as a float 174500.12. This code is ugly but parses hms or dms:

from astropy.coordinates import Angle

def to_sigproc_angle(angle_val):
    """ Convert an astropy.Angle to the ridiculous sigproc angle format string. """
    x = str(angle_val)

    if '.' in x:
        if 'h' in x:
                d, m, s, ss = int(x[0:x.index('h')]), int(x[x.index('h')+1:x.index('m')]), \
                int(x[x.index('m')+1:x.index('.')]), float(x[x.index('.'):x.index('s')])
        if 'd' in x:
            d, m, s, ss = int(x[0:x.index('d')]), int(x[x.index('d')+1:x.index('m')]), \
            int(x[x.index('m')+1:x.index('.')]), float(x[x.index('.'):x.index('s')])
    else:
        if 'h' in x:
            d, m, s = int(x[0:x.index('h')]), int(x[x.index('h')+1:x.index('m')]), \
            int(x[x.index('m')+1:x.index('s')])
        if 'd' in x:
            d, m, s = int(x[0:x.index('d')]), int(x[x.index('d')+1:x.index('m')]), \
            int(x[x.index('m')+1:x.index('s')])
        ss = 0
    num = str(d).zfill(2) + str(m).zfill(2) + str(s).zfill(2)+ '.' + str(ss).split(".")[-1]
    return np.float64(num)

...

# Get RA and DEC and convert to sigproc format
ra = Angle(ihdr['raj'], unit='hourangle')
dec = Angle(ihdr['dej'], unit='degree')
sigproc_hdr['src_raj'] =  to_sigproc_angle(ra)
sigproc_hdr['src_dej'] =  to_sigproc_angle(dec)

Crash bug

The code fails on L204 and L291 as axnames is of type list, while the RHS a tuple.

Easiest fix is:

 if ndim >= 3 and axnames[-3:] == ('time', 'pol', 'freq'):
becomes
 if ndim >= 3 and tuple(axnames[-3:]) == ('time', 'pol', 'freq'):
@benbarsdell
Copy link
Collaborator

Thanks for finding these bugs.

Couldn't you use http://docs.astropy.org/en/stable/api/astropy.coordinates.Angle.html#astropy.coordinates.Angle.dms

Astropy does seem like a heavyweight solution for this. How/where do you need to interact with the angle value in the headers? E.g., are you setting them in a source block, loading them from a Sigproc read block, something else? I guess we need to decide on a standard representation of them in headers.

@telegraphic
Copy link
Collaborator Author

Yes, generally read from guppi raw files, which are often strings (unsure what the spec is, or if people actually follow it).

An yes, that code is particularly horrible (not written by me in my defense!), and came from a non-astropy-using codebase.

Astropy angles and units (http://docs.astropy.org/en/v0.2.1/units/index.html) are both pretty useful. I think we either use them extensively (instead of pint?), or not at all...

@jaycedowell
Copy link
Collaborator

@telegraphic is this still a thing? It looks like the crash has been fixed but nothing was done about the angle representation.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Development

No branches or pull requests

3 participants