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

add make_pretty_from_fits file to develop #319

Merged
merged 23 commits into from
Jan 16, 2018
Merged

add make_pretty_from_fits file to develop #319

merged 23 commits into from
Jan 16, 2018

Conversation

jermainegug
Copy link
Contributor

@jermainegug jermainegug commented Jan 12, 2018

*creates norm variable
*checks wcs (**issue)
*changed cmap from cubealpha_r to inferno

try: wcs = WCS(fname) except: wcs = False
This does not work as expected.
If we do:

wcs = WCS(fname) print(wcs)

We get this as a result:

Number of WCS axes: 2
CTYPE : ''  ''  
CRVAL : 0.0  0.0  
CRPIX : 0.0  0.0  
PC1_1 PC1_2  : 1.0  0.0  
PC2_1 PC2_2  : 0.0  1.0  
CDELT : 1.0  1.0  
NAXIS : 3352  2532

*creates norm variable
*checks wcs (issue)
*changed cmap from cubealpha_r to inferno
@wtgee
Copy link
Member

wtgee commented Jan 12, 2018

Thanks @jermainegug. So when I told you to develop it I had you work with a separate script because you didn't have POCS installed yet. But know that you do...

What we want to do is take your _make_pretty_from_fits and repalce what is in: https://github.com/panoptes/POCS/blob/develop/pocs/utils/images/__init__.py#L87. Then you can just go ahead and delete this file.

@wtgee
Copy link
Member

wtgee commented Jan 12, 2018

You can just keep working in the same branch and any new changes you push will automatically become part of this PR.

*using the function in make_pretty_from_fits file
*delete make_pretty_from_fits file
@AnthonyHorton
Copy link
Collaborator

Ah, OK, going to need to try something else to check for a valid WCS. Seems WCS(fname) will just create a mostly empty template in the absence of a full WCS in the FITS headers. Check the astropy.wcs documentation, there must be a function for checking whether a WCS is present/valid or not.

Also, for the plot using WCS use RA & dec for the axes (not an overlay) and don't bother with Galactic coordinates, at least not by default. Could make the choice of coordinate systems an option, of course.

@wtgee
Copy link
Member

wtgee commented Jan 12, 2018

Just about to put the comment in review for WCS. We do this elsewhere. Basically:

w = wcs.WCS(filename)
assert w.is_celestial

wcs = False

"""
try:
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

try:
    w = wcs.WCS(filename)
    assert w.is_celestial
except AssertionError:
    pass

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Or just skip the try and use if wcs.is_celstial below.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Get this error message AttributeError: 'bool' object has no attribute 'WCS' for w = wcs.WCS(filename).

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@jermainegug That's because in your code you're setting wcs to False, whereas @wtgee's code suggestion is assuming that wcs is astropy.wcs (i.e. you're doing from astropy import wcs)

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes, I didn't adjust it for your code, I just copied what we have elsewhere.

You import WCS directly so use that instead of wcs.WCS.

But, you are also mixing types. You create wcs = False (a boolean) up above and then assign wcs=WCS(filename) (a WCS class) to it. I would just skip the try and use the is_celestial in the if statement directly after getting the wcs.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

What is wcs in this case? wcs=WCS(fname)?

Copy link
Collaborator

@AnthonyHorton AnthonyHorton Jan 12, 2018

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

R.T.D.

It's a WCS object, i.e. an instance of the astropy.wcs.WCS class.

http://docs.astropy.org/en/stable/api/astropy.wcs.WCS.html#astropy.wcs.WCS

In retrospect, probably best to use a more distinct variable name for it, e.g. image_wcs.

@@ -6,7 +6,10 @@
from matplotlib import pyplot as plt
from warnings import warn

from astropy.io import fits
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I don't think you are using this.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Nope, using that in my jupyter notebook.

*fixed check for wcs
*not sure about RA & dec axes
@codecov
Copy link

codecov bot commented Jan 12, 2018

Codecov Report

Merging #319 into develop will decrease coverage by 0.75%.
The diff coverage is 83.87%.

Impacted file tree graph

@@             Coverage Diff             @@
##           develop     #319      +/-   ##
===========================================
- Coverage    83.81%   83.06%   -0.76%     
===========================================
  Files           51       51              
  Lines         3367     3454      +87     
  Branches       422      439      +17     
===========================================
+ Hits          2822     2869      +47     
- Misses         404      440      +36     
- Partials       141      145       +4
Impacted Files Coverage Δ
pocs/utils/images/__init__.py 82.75% <83.87%> (-2.68%) ⬇️
pocs/state/states/default/tracking.py 83.33% <0%> (-16.67%) ⬇️
pocs/utils/rs232.py 79.64% <0%> (-13.98%) ⬇️
pocs/scheduler/observation.py 93.44% <0%> (-6.56%) ⬇️
pocs/dome/astrohaven.py 72.38% <0%> (-6.29%) ⬇️
pocs/utils/database.py 89.13% <0%> (-4.35%) ⬇️
pocs/observatory.py 84.84% <0%> (-1.36%) ⬇️
pocs/state/states/default/pointing.py 87.75% <0%> (ø) ⬆️
pocs/scheduler/field.py 100% <0%> (ø) ⬆️
pocs/mount/simulator.py 92.79% <0%> (+0.26%) ⬆️
... and 2 more

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update f90359a...3c7ebf9. Read the comment docs.

*just set axis labels to RA and dec
@AnthonyHorton
Copy link
Collaborator

AnthonyHorton commented Jan 12, 2018

@jermainegug Have a look through the documentation on this under the astropy.visualization module: http://docs.astropy.org/en/stable/visualization/wcsaxes/index.html

For this you should be using matplotlibs object orientated interface, rather than the pyplot interface, see the second example in the above link.

That will then make it easier to properly format the tick labels for axis, see: http://docs.astropy.org/en/stable/visualization/wcsaxes/ticks_labels_grid.html You're going to want to do this in order to make the Right Ascension axis use hour units instead of the default degrees.

Also, you need to make the axis labels more explicit. I'd prefer to spell out Right Ascension & declination (if there's room) but more importantly you should specify the coordinate system and (if applicable) equinox that are being used. The WCS object should have this information (and also allow you to transform to others). It's probably either FK5 J2000 or ICRS.

EDIT: There's an example of setting tick label format via the object orientated interface in the Huntsman dither.py utils: https://github.com/AstroHuntsman/huntsman-pocs/blob/develop/huntsman/utils/dither.py

*descriptive axes
*coordinate system in title?
*how to get equinox?
@@ -89,7 +90,7 @@ def make_pretty_image(fname, timeout=15, **kwargs): # pragma: no cover
def _make_pretty_from_fits(fname, **kwargs):
config = load_config()

title = '{} {}'.format(kwargs.get('title', ''), current_time().isot)
title = '{} {} {}'.format(kwargs.get('title', ''), current_time().isot, 'WCS') # is WCS the coordinate system?
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

You can do current_time(pretty=True) for a better time output here.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

And WCS = "world coordinate system", but you don't just want the string there unless you are only putting it there when is_celestial is true or something.

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Even then it will be obvious from the rest of the plot whether a WCS from the FITS image has been used or not, from the axis labels, tick labels and grid. Don't need to put it in the title too.

Copy link
Contributor Author

@jermainegug jermainegug Jan 12, 2018

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

So we will just be able to tell the coordinate system based off what the image looks like? Thus no need to have anything writing out that it uses WCS.

ax.set_xlabel('RA')
ax.set_ylabel('Dec')

ra = ax.coords[0]
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This is not really an ra, it is an ra_axis. Same with the dec. Better to be explicit about the name because if you had another variable (which maybe you will in future) that actually holds an RA value you would get confused by this.

ra = ax.coords[0]
dec = ax.coords[1]

ra.set_axislabel('Right Ascension')
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I tend to do label [units] (although I think @AnthonyHorton does something different), in which case this would be Right Ascension [hms] and Declination [dms].

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes, I tend to prefer the label / units convention (i.e. the numbers are the quantity in question divided by the given unit). Either is fine, though.

In this particular case putting units in the axis labels should be redundant though, the tick label formats include units, see https://github.com/AstroHuntsman/HuntsmanPOCSLegacy/blob/develop/examples/notebooks/dice9.png for an example.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I have written them in for now. Will make sure that the units are in the axis label when I run the code with a fits file that has wcs., then I will get rid of the units in the label.

ax.set_xlabel('Pixel x-axis')
ax.set_ylabel('Pixel y-axis')

ax.coords[0].set_axislabel('Pixel x-axis')
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

As per comments above, I would just do X [pixels] and Y [pixels].

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Or X / pixels, Y / pixels ;-)

@@ -103,14 +104,22 @@ def _make_pretty_from_fits(fname, **kwargs):
ax = plt.subplot(projection=wcs)

ax.coords.grid(True, color='white', ls='solid')
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This line could be moved below if clause so it is not repeated.

dec.set_axislabel('Declination')

ra.set_major_formatter('hh:mm:ss')
dec.set_major_formatter('dd:mm:ss')
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

For Huntsman, PANOPTES or in fact most other optical instruments the field of view will be large enough that including seconds in the major tick labels is unnecessary. hh:mm and dd:mm is probably better.


new_filename = fname.replace('.fits', '.jpg')

data = getdata(fname)
plt.imshow(data, cmap='cubehelix_r', origin='lower')

norm = ImageNormalize(interval=MinMaxInterval(), stretch=LogStretch())
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Might be better to set the interval to trim off a (very) small percentage at either end, so that the normalisation doesn't get driven by dead pixels, hot pixels, cosmic rays and saturated stars.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

What sort of percentage are we thinking? Something like 0.1%?

*for lines with (?) contents should be added to config?
*change from minmax to percentile interval
*add alpha option for grid (transparent lines)
wcs = WCS(fname)
ax = plt.subplot(projection=wcs)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This line might need to stay inside the if. I'm not sure what happens if you pass a projection but don't have a valid WCS. Have you tried on the original images that weren't solved yet?

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Just for fun I just ran it on an image that a PANOPTES unit just took that has a pretty empty star field and a meteor streak. Looks like it handles the lack of wcs fine.

test

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

There's a lot of dead white space in there around the edges. I think you can do a fig.tight_layout() to get rid of some of that.

Copy link
Contributor Author

@jermainegug jermainegug Jan 12, 2018

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I tested it on my Jupyter notebook and it seems to work just fine. If I don't assign anything to projection then I can't use coords[] for the pixel coordinates. I would have to use this for else
ax = plt.subplot()
ax.set_xlabel('X / pixels')
ax.set_xlabel('Y / pixels')

But I think what I have should work just fine :)
If not then we can change it to this.

Awesome, image looks good, is this with percent_value = 99.81? Only asking because maybe we should make that value able to be modified since when I was trying different percentages, and even changing it by 0.1% it changes quite drastically.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

And I will add the fig.tight_layout() soon and let you know how it looks.

@wtgee
Copy link
Member

wtgee commented Jan 12, 2018

This has nothing to do with your code as it was already in the function, but we shouldn't actually be doing the symlink stuff here since we don't do it in _make_pretty_from_cr2. Can you remove all those lines since you are updating the function anyway? Thanks!

@@ -87,12 +90,34 @@ def make_pretty_image(fname, timeout=15, **kwargs): # pragma: no cover
def _make_pretty_from_fits(fname, **kwargs):
config = load_config()

title = '{} {}'.format(kwargs.get('title', ''), current_time().isot)
title = '{} {}'.format(kwargs.get('title', ''), current_time(pretty=True).isot)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Tests are failing on this. Doing the pretty=True turns it into a string, so then you can't call .isot on it. You want one or the other. They are basically the same except pretty has spaces for easier reading.

@wtgee
Copy link
Member

wtgee commented Jan 12, 2018 via email

@jermainegug
Copy link
Contributor Author

horsehead_inferno_new_method

I used the example image and the wcs is working great! Like @AnthonyHorton mentioned, we wont need the units in the axis label since they are written to values on the image.

@wtgee
Copy link
Member

wtgee commented Jan 13, 2018

Note that you are failing on some codestyle issues. See https://travis-ci.org/panoptes/POCS/builds/328330892#L4924

With atom or sublime you can use an autopep8 tool that should correct most things.

@jermainegug
Copy link
Contributor Author

autopep8 doesn't change anything in my file..

@wtgee
Copy link
Member

wtgee commented Jan 13, 2018

autopep8 doesn't change anything in my file..

They were related to the comments you deleted anyway.

Can you post a pic with the tight_layout? Also, in pick above I see the title but not a date. Seems like there should be both from this code.

@jermainegug
Copy link
Contributor Author

jermainegug commented Jan 15, 2018

Check seems to give us this:

pocs/tests/utils/test_image_utils.py::test_clean_observation_dir
  /home/travis/build/panoptes/POCS/pocs/utils/images/__init__.py:140: 
UserWarning: Can't link latest image: [Errno 2] No such file or directory: 
'/tmp/tmppbxfbjvv/solved.jpg' -> '/var/panoptes/images/latest.jpg'
    warn("Can't link latest image: {}".format(e))

@wtgee
Copy link
Member

wtgee commented Jan 15, 2018

See my comment above about removing all the link stuff: #319 (comment)

@jermainegug
Copy link
Contributor Author

Oh yes, I forgot to delete those. Thanks.

@wtgee
Copy link
Member

wtgee commented Jan 15, 2018

Oh yes, I forgot to delete those. Thanks.

There is other stuff that goes with it. For instance once you delete that I don't think you need the config or the image_dir. Make sure to do a full cleanup please. Thanks!

@wtgee
Copy link
Member

wtgee commented Jan 15, 2018

Regarding the title and time from the fits headers, yes I would think so. The example I showed up lets a title in the kwargs override that so it would just be a backup.

@jermainegug
Copy link
Contributor Author

Ignore the recent commit.

@wtgee
Copy link
Member

wtgee commented Jan 15, 2018

You want to keep the kwargs. That lets someone pass something in to override. So if I am using make_pretty_image interactively from a jupyter console (which I do freqently) I can say make_pretty_image(title="test_01") and have it work. That should take precedence over the fits header because I am explicitly asking for it.

@wtgee
Copy link
Member

wtgee commented Jan 15, 2018

Ignore the recent commit.

Too late. :)

(?) I am not sure if I should use the imported current time or use the DATE-OBS from the header.
*change title

image_dir = config['directories']['images']
percent_value = 99.9
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Just to make it flexible for the future I would pull this from kwargs too. So kwargs.get('normalize_clip_percent', 99.9)

ax.set_xlabel('X / pixels')
ax.set_ylabel('Y / pixels')

ax.imshow(data, norm=norm, cmap='inferno', origin='lower')
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Also make this slightly more flexible by doing a cmap=kwargs.get('cmap', 'inferno'). Then you can play with different colormaps in runtime without touching code.

*allow for selection of percent, color map
*neaten up code
@jermainegug
Copy link
Contributor Author

Yay! All checks passed :)

Copy link
Member

@wtgee wtgee left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

One more small change and then I think it's good.


title = '{} {}'.format(kwargs.get('title', ''), current_time().isot)
title = kwargs.get('title', header.get('FIELD', ''))
date_time = header.get('DATE-OBS', current_time()).replace('T', ' ', 1)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This might need to be tested. current_time will return a time object and then the replace will probably fail. If you leave at current_time(pretty=True) it will return a string and the replace will work on a string (it just won't replace anything since the T has already been removed, but that's fine), which is probably safer and means you don't need to bother testing that.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I have changed it now but I will test anyways (want to make sure that it doesn't get rid of the T in UTC).

Copy link
Contributor Author

@jermainegug jermainegug Jan 15, 2018

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Oh, current_time doesn't print out UTC. So it is all good then. Works just fine.

Copy link
Member

@wtgee wtgee left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

LGTM! We've beat this PR to death so any additional items can be fixed via an Issue.

@jermainegug I'll let you do the honors of merging your own PR. Make sure to select the Squash option and to say "Fixes #232" in merge commit. You are also okay to delete the branch after (it will show you a button). Thanks!

@jermainegug
Copy link
Contributor Author

jermainegug commented Jan 16, 2018

Indeed we have!

No problem but I do not have write access so the merge button does not show up.

Not sure if you can do this but did you want me to say it also fixes #43 along with #232 ?
And I guess we can close AstroHuntsman#2

@wtgee
Copy link
Member

wtgee commented Jan 16, 2018

Ok, I'll go ahead and merge. #43 stays open as it is more of a question of what we want to do with the pretty pictures than it is about getting pretty pictures working.

@wtgee wtgee merged commit 46d7ae8 into panoptes:develop Jan 16, 2018
@jermainegug jermainegug deleted the make-pretty-from-fits branch January 16, 2018 22:09
@wtgee wtgee mentioned this pull request Jan 24, 2018
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

3 participants