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 function to plot Fermi catalog light curves #286

Merged
merged 29 commits into from Jun 29, 2015

Conversation

@JonathanDHarris
Copy link
Contributor

@JonathanDHarris JonathanDHarris commented Jun 13, 2015

@cdeil
Here's a first version for review. I'm not that familiar with Python, so don't be afraid to be pedantic about:
-naming conventions,
-how / where I import things
-comment styles
-etc

At the moment the tool plots Fermi MET as the x-axis. I don't know if there is a package default to use MJD or something else, in which case I can change it. I'm also just showing the plot when you call the function, if you want an object returned then let me know, I might need an example of the format.

@cdeil cdeil added the feature label Jun 13, 2015
@cdeil cdeil added this to the 0.3 milestone Jun 13, 2015
@cdeil cdeil self-assigned this Jun 13, 2015
@@ -2,31 +2,126 @@
from __future__ import (absolute_import, division, print_function,
unicode_literals)

__all__ = ['plot_time_difference_distribution',
from ..datasets import fetch_fermi_catalog

This comment has been minimized.

@cdeil

cdeil Jun 13, 2015
Member

This fails on travis-ci:
https://travis-ci.org/gammapy/gammapy/jobs/66665216#L1200
It should fail locally for you too if you run

python setup.py test -P time

I think the issue is circular imports (see #171) and the short-term solution is to move the import into the plot_light_curve function.

]


def plot_time_difference_distribution(time, ax=None):

This comment has been minimized.

@cdeil

cdeil Jun 13, 2015
Member

The diff on Github looks complex because you've removed the plot_time_difference_distribution function.
Can you please put that back in, I wanted to implement it eventually.

The easiest might be to copy & paste it from here:
https://github.com/gammapy/gammapy/blob/master/gammapy/time/plot.py#L9

Axes
name_3FGL : `string`
The 3FGL catalog name of the object to plot
fermi_met_start : `int`

This comment has been minimized.

@cdeil

cdeil Jun 13, 2015
Member

I don't know if you've seen it, yesterday I wrote some notes on time handling in Gammapy

Any of these would be valid choices I think:

  • Don't take start and end time as input, always plot the full light curve.
  • Take two astropy.time.Time objects as input called start_time and end_time
  • Take one astropy.time.Time object as input with length 2 called time_range
End time of the light curve in Fermi MET
Usage
gammapy.time.plot_light_curve('3FGL J0349.9-2102', 2.1e8, 3.2e8)

This comment has been minimized.

@cdeil

cdeil Jun 13, 2015
Member

You should add a full example to the docstring.
This means adding a section:

Examples
-------------

and a plot directive as shown e.g. in this example: https://gammapy.readthedocs.org/en/latest/api/gammapy.irf.EffectiveAreaTable.html

So Sphinx will actually execute that code and put the image in the docs.
Fetching lots of datasets in the docs build will become problematic, but for now I think it's OK.


try:
catalog_index = np.where(fermi_cat[1].data['Source_Name'] == name_3FGL)[0][0]
except IndexError:

This comment has been minimized.

@cdeil

cdeil Jun 13, 2015
Member

You should remove all of the try ... except in this function.

At least in library code like you're writing here it's important to raise Exceptions or Warnings when something goes wrong. Just printing an error message means that callers don't have a chance to catch the issue if they e.g. produce 1000s or plots in a pipeline and one source name is wrong.

Of course sometimes it can make sense to catch errors and re-raise more specific errors, or to catch warnings and silence them completely because it's expected behaviour in that case, but it's rare.

Plus we can't write Gammapy so that we check inputs all the time.
It would literally mean writing ~ 2x as much code.
See some comments I made here: https://gammapy.readthedocs.org/en/latest/development/index.html#what-checks-and-conversions-should-i-do-for-inputs


upper_lims_y = np.copy(flux_history)

upper_lims_x = \

This comment has been minimized.

@cdeil

cdeil Jun 13, 2015
Member

In Gammapy (and most other Python projects these days I think), please don't use \.
Here you can simply put the statement on one line.

If you have long expressions, you can use parenthesis:

result = (aaa *
             bbb)

But I much prefer the pattern where one introduces temp variables.
In this case

idx = np.where(flux_history_lower_bound <= 0)

would be a good temp variable because you re-use it several times, so you actually avoid duplicated code.

# TODO: implement!
raise NotImplementedError
# Plot data points and upper limits.
plt.errorbar(time_mid, flux_history, yerr=(flux_history_lower_bound, flux_history_upper_bound), \

This comment has been minimized.

@cdeil

cdeil Jun 13, 2015
Member

No need for \ at the end of the lines. The fact that the plt.errorbar() function call has a start and stop parenthesis means it's OK to use multiple lines for the function arguments.

plt.errorbar(upper_lims_x[np.where(upper_lims_y <= 0)], upper_lims_y[np.where(upper_lims_y <= 0)], \
yerr=(upper_lims_y_end[np.where(upper_lims_y <= 0)], upper_lims_y_start[np.where(upper_lims_y <= 0)]), \
marker=None, elinewidth=1, linewidth=0, lolims=True, color='black')
plt.xlabel('time [Fermi MET]')

This comment has been minimized.

@cdeil

cdeil Jun 13, 2015
Member

I would much prefer to get a time axis with years and maybe months as labels than Fermi MET.
Those numbers mean nothing to me and probably to most other astronomers as well.

Here's some info how to convert times to a format matplotlib understands:
http://astropy.readthedocs.org/en/latest/time/index.html?highlight=matplotlib#get-representation

This will be the hard part of this PR ... figuring out how to correctly create astropy.timeTime objects and how to correctly plot them with matplotlib.
If you have any questions or get stuck, please ask and I can investigate how to do this also.

@cdeil
Copy link
Member

@cdeil cdeil commented Jun 13, 2015

@JonathanDHarris - I've left inline comments. Please ask if it's unclear what I mean or what should be done. I'll probably have some more style comments later, so please be patient ... review of the first few pull requests for new contributors usually takes a while with any project. For now the important thing is to have a working example in the docstring and, if you want to take this on, have a correct time axis with years and months.

@cdeil cdeil changed the title First version of a Fermi light curve tool Add function to plot Fermi catalog light curves Jun 13, 2015
@cdeil
Copy link
Member

@cdeil cdeil commented Jun 13, 2015

One more thing: you can drag & drop PNG images into the Github comment field.
Could you please do this for this PR?

@cdeil
Copy link
Member

@cdeil cdeil commented Jun 13, 2015

@JonathanDHarris - I think plotting functions should take an axes object as input and also return it.
See astropy/astroplan#7 (comment)
Once the discussion there is settled, we'll document this pattern here.

It's quite a bit of boilerplate code, but it means a lot of flexibility for the caller, e.g. users could easily plot several light curves in panels (example).

@JonathanDHarris
Copy link
Contributor Author

@JonathanDHarris JonathanDHarris commented Jun 14, 2015

@cdeil
I think I've taken into account all of your comments. I do need to go back through and double check. I'm on holiday this week so I might be out of touch for a little while.

python setup.py test -P time shows 4 tests passed for me now.

I had trouble pulling from origin master earlier, but I don't think it has effected this work.

@JonathanDHarris
Copy link
Contributor Author

@JonathanDHarris JonathanDHarris commented Jun 14, 2015

exampleplot


def plot_time_difference_distribution(time, ax=None):
"""Plot event time difference distribution.

This comment has been minimized.

@cdeil

cdeil Jun 14, 2015
Member

Put the empty line back in.
(I think it's actually required for the docstring to be processed correctly by Sphinx)

@@ -30,3 +28,140 @@ def plot_time_difference_distribution(time, ax=None):

# TODO: implement!
raise NotImplementedError

def plot_light_curve(name_3FGL, time_start, time_end, ax = None):

This comment has been minimized.

@cdeil

cdeil Jun 14, 2015
Member

Rename function to the more explicit plot_fermi_3fgl_light_curve.

Add default values time_start=None and time_end=None and document / implement them to mean "full available time range".

This comment has been minimized.

@cdeil

cdeil Jun 14, 2015
Member

Like most Python project we follow the PEP8 coding style.
The pep8 tool tells you what you need to fix, on this line you should write ax=None instead of ax = None, i.e. remove the spaces around the equals sign:

$ pep8 gammapy/time/plot.py 
gammapy/time/plot.py:9:1: E302 expected 2 blank lines, found 1
gammapy/time/plot.py:32:1: E302 expected 2 blank lines, found 1
gammapy/time/plot.py:32:57: E251 unexpected spaces around keyword / parameter equals
gammapy/time/plot.py:32:59: E251 unexpected spaces around keyword / parameter equals
gammapy/time/plot.py:63:80: E501 line too long (91 > 79 characters)
gammapy/time/plot.py:81:80: E501 line too long (92 > 79 characters)
gammapy/time/plot.py:84:80: E501 line too long (81 > 79 characters)
gammapy/time/plot.py:86:80: E501 line too long (89 > 79 characters)
gammapy/time/plot.py:89:80: E501 line too long (90 > 79 characters)
gammapy/time/plot.py:91:80: E501 line too long (82 > 79 characters)
gammapy/time/plot.py:102:80: E501 line too long (101 > 79 characters)
gammapy/time/plot.py:104:80: E501 line too long (120 > 79 characters)
gammapy/time/plot.py:105:80: E501 line too long (120 > 79 characters)
gammapy/time/plot.py:146:80: E501 line too long (90 > 79 characters)
gammapy/time/plot.py:147:80: E501 line too long (98 > 79 characters)
gammapy/time/plot.py:150:80: E501 line too long (99 > 79 characters)
gammapy/time/plot.py:151:18: E265 block comment should start with '# '
gammapy/time/plot.py:153:80: E501 line too long (100 > 79 characters)
gammapy/time/plot.py:154:80: E501 line too long (117 > 79 characters)
gammapy/time/plot.py:155:80: E501 line too long (83 > 79 characters)
gammapy/time/plot.py:156:80: E501 line too long (102 > 79 characters)
gammapy/time/plot.py:157:80: E501 line too long (119 > 79 characters)
gammapy/time/plot.py:158:80: E501 line too long (84 > 79 characters)
gammapy/time/plot.py:165:5: E265 block comment should start with '# '
gammapy/time/plot.py:167:14: W292 no newline at end of file

You can install and run the autopep8 tool to automatically fix those.
One of the few exceptions is that sometimes it makes sense to have lines > 80 chars, so my suggestion would be to always run git diff before accepting all autopep8 suggestions.

You could also give this Python IDE a try: https://www.jetbrains.com/pycharm/
Open-source projects like Gammapy get free licenses for team members, so if you want one, please email me.

name_3FGL : `string`
The 3FGL catalog name of the object to plot
fermi_met_start : `~astropy.time.Time`
Astropy time object for the start of the light curve

This comment has been minimized.

@cdeil

cdeil Jun 14, 2015
Member

No need to write "Astropy time object" again, that's what the ~astropy.time.Time in the previous line says.
Here I'd write:

start_time : `~astropy.time.Time`
    Light curve start time

i.e. also change the parameter name, because the user no longer has to give MET.

@JonathanDHarris
Copy link
Contributor Author

@JonathanDHarris JonathanDHarris commented Jun 19, 2015

exampleplot

@cdeil
Copy link
Member

@cdeil cdeil commented Jun 19, 2015

There's still this issue on travis-ci:
https://travis-ci.org/gammapy/gammapy/jobs/67544898#L1856

@JonathanDHarris – Are you seeing the same locally?

I'll try it out locally now.

@cdeil
Copy link
Member

@cdeil cdeil commented Jun 19, 2015

I don't know why I didn't think about this before, but now I think it would be better to refactor this code into a LightCurve class that sub-classes Table and has a plot method. Then this could be re-used for other light curves, e.g. from HESS or X-ray or whatever.
However I'm not suggesting you try to do this here ... creating the LightCurve class will take thought and discussion and more work (e.g. for FITS / CSV I/O) ... just wanted to mention it here because I thought of it just now that this will have to be done in the future.

+ astropy.time.TimeDelta(upper_lims_x, format='sec')).plot_date

# Plot data points and upper limits.
plt.errorbar(time_mid, flux_history,

This comment has been minimized.

@cdeil

cdeil Jun 22, 2015
Member

You should explicitly draw on ax, i.e. say ax.errorbar, ax.xlabel, ...

This comment has been minimized.

@JonathanDHarris

JonathanDHarris Jun 27, 2015
Author Contributor

Sorry, I don't understand. I tried this but I got
AttributeError: 'AxesSubplot' object has no attribute 'xlabel'

This comment has been minimized.

@cdeil

cdeil Jun 27, 2015
Member

Sorry ... it's not ax.xlabel ... it's ax.set_xlabel.

marker=None, elinewidth=1, linewidth=0, lolims=True, color='black')
plt.xlabel('date')
plt.ylabel('flux [ph/cm^2/s]')
plt.gca().xaxis_date()

This comment has been minimized.

@cdeil

cdeil Jun 22, 2015
Member

Same here ... use ax instead of plt.gca().

time_start = Time('2010-01-01T00:00:00')
time_end = Time('2015-02-02T02:02:02')
plt.plot = plot_fermi_3fgl_light_curve('3FGL J0349.9-2102', time_start, time_end)

This comment has been minimized.

@cdeil

cdeil Jun 22, 2015
Member

For me this example looks like this:

image

The display of upper limits with the arrows pointing up is confusing, no?

They've probably discussed this at length ... maybe it's best to just plot it the way it's done in the 3FGL paper?
screen shot 2015-06-22 at 09 59 48

This comment has been minimized.

@cdeil

cdeil Jun 22, 2015
Member

(except for the time axis .... for this I still like years, months a lot, so please keep this, at least as an option)

This comment has been minimized.

@cdeil

cdeil Jun 27, 2015
Member

You shouldn't assign to plt.plot ... that's a function in matplotlib that you're overwriting with the axes object returned here.
Here there's no need to store the axes object in a variable, because you don't want to modify the plot further, so just remove plt.plot =.

This comment has been minimized.

@JonathanDHarris

JonathanDHarris Jun 27, 2015
Author Contributor

This may also have been causing the Travis issue, we'll see.

@cdeil
Copy link
Member

@cdeil cdeil commented Jun 22, 2015

I can't reproduce the Sphinx errors TypeError: 'AxesSubplot' object is not callable locally:
https://travis-ci.org/gammapy/gammapy/jobs/67544898#L1855
It's unrelated to this PR ... if it doesn't go away with the next travis-ci run here I'll investigate.

This PR is almost ready.

One thing that would be nice is if you squash the 24 commits here down into a single commit "Add function to plot Fermi 3FGL light curves". I've usually used git rebase -i master as described here, but in this case it might also make sense to do this non-interactively as described here.
If you decide to do this, please make a backup copy of your repo before you do. If you don't want to, that's fine, too, a clean git commit history of Gammapy is not terribly important.

@JonathanDHarris
Copy link
Contributor Author

@JonathanDHarris JonathanDHarris commented Jun 26, 2015

@cdeil
It's very odd that the arrow direction on the limit is the other way around for you. This is what I see. Any idea why that might be?

exampleplot

@JonathanDHarris
Copy link
Contributor Author

@JonathanDHarris JonathanDHarris commented Jun 26, 2015

I'm changing the style to match the one you suggested, because the pyplot upperlimits are always fiddly. But if you know why we saw them point in opposite directions I'd be interested to know for future reference.

exampleplot

plt.show()
"""
from ..datasets import fetch_fermi_catalog
from ..time.utils import TIME_REF_FERMI

This comment has been minimized.

@cdeil

cdeil Jun 27, 2015
Member

Can you move the TIME_REF_FERMI import to the top of the file?
That should be the default, and only if that doesn't work for some reason (circular import, optional dependency), then it's moved into functions / methods.

@cdeil
Copy link
Member

@cdeil cdeil commented Jun 27, 2015

@JonathanDHarris – Now there's a merge conflict (probably in the changelog and easy to resolve), which is the reason travis-ci isn't running. Can you please rebase on master?

@JonathanDHarris
Copy link
Contributor Author

@JonathanDHarris JonathanDHarris commented Jun 28, 2015

Hm. Travis fails again. This is causing hard to trouble shoot because if I run python setup.py build_sphinx the build fails, but if I run exactly the same command again it will then pass.

@cdeil
Copy link
Member

@cdeil cdeil commented Jun 28, 2015

@JonathanDHarris – The reason for the travis-ci fails is: astropy/astroplan#15

Could you please uncomment this line?
https://github.com/gammapy/gammapy/blob/master/gammapy/conftest.py#L35
Then they should pass and we can merge this.

For the future ... you should have rebased on master, not merged master into your branch:
https://github.com/gammapy/gammapy/network
That way Gammapy will have a simple, linear history.
For this time, it doesn't matter ... leave as is.

ax.xaxis_date()
ax.xaxis.set_major_formatter(mdates.DateFormatter('%m/%Y'))
ax.xaxis.set_major_locator(mdates.MonthLocator(interval=6))
plt.gcf().autofmt_xdate()

This comment has been minimized.

@cdeil

cdeil Jun 28, 2015
Member

To stay fully object-oriented and avoid side effects, you should access the correct figure via the axes:

ax.figure.autofmt_xdate()

This comment has been minimized.

@JonathanDHarris

JonathanDHarris Jun 28, 2015
Author Contributor

Thanks, I was looking for a while how to do that.

@cdeil
Copy link
Member

@cdeil cdeil commented Jun 28, 2015

This now looks good to me and locally this works for me also.
Let's hope travis-ci passes, then I'll merge tonight.

@JonathanDHarris
Copy link
Contributor Author

@JonathanDHarris JonathanDHarris commented Jun 28, 2015

Fingers crossed!

@cdeil
Copy link
Member

@cdeil cdeil commented Jun 29, 2015

@JonathanDHarris – Another fail came up that's unrelated to this PR:
https://travis-ci.org/gammapy/gammapy/jobs/68707705#L1115
It's a recently-introduced bug in Astropy:
astropy/astropy#3891 (comment)

I'll disable that test for now in Gammapy with a commit directly to master and then re-start the build here. Then it should pass.

@cdeil
Copy link
Member

@cdeil cdeil commented Jun 29, 2015

The travis-ci build for master is broken here anyways because of this issue that astroplan downloads a file on import:
https://travis-ci.org/gammapy/gammapy/jobs/68773866#L1135

So I'm merging this now, and if there's some other issue I'll fix it in master.
I don't expect there to be though ... this looks good!

Thanks for the nice contribution and for being patient with code review and dealing with unrelated issues!

cdeil added a commit that referenced this pull request Jun 29, 2015
Add function to plot Fermi catalog light curves
@cdeil cdeil merged commit 6553d24 into gammapy:master Jun 29, 2015
1 check was pending
1 check was pending
continuous-integration/travis-ci/pr The Travis CI build is in progress
Details
@cdeil
Copy link
Member

@cdeil cdeil commented Jun 29, 2015

@JonathanDHarris – I've decided to offer the user the possibility to optionally pass in strings to specify the time to make this simpler ... I just forward those along to the Time constructor: 9c979c5

And I've opened #294 as a TODO for me that we need to set up testing of Figures.

@JonathanDHarris
Copy link
Contributor Author

@JonathanDHarris JonathanDHarris commented Jun 29, 2015

@cdeil Awesome! Thanks for all your help on this.

@JonathanDHarris JonathanDHarris deleted the JonathanDHarris:FermiLightCurves branch Jul 30, 2016
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Linked issues

Successfully merging this pull request may close these issues.

None yet

2 participants