Extending mainpage example - #114

This is what the image looks like at the moment:

The dark green is not well visible ... maybe just show all apertures in the light red?
And maybe increase the circle radius by a pixel or two to make the apertures around the bright sources better visible?

@cdeil - what about using an inverted grayscale for the image?

Inverted grayscale could look nicer ... @bsipocz can you try and paste the result here before we decide?

Yes, I guess the inverted gray looks way much better.

@bsipocz - cool, can you upload an example here? (I mean in a comment, so I can see)

@astrofrog - yes, I'm trying, but can't figure out how to do it via github as @cdeil 's figure.

Anyway here it is:
http://star.herts.ac.uk/~bsipocz/astropy/index-1.png

@bsipocz You can drag & drop images into the Github issue input fields: https://github.com/blog/1347-issue-attachments

Wow, thanks @cdeil.

I like the inverted color scale.
This is going to be used e.g. in presentations as an example image, so I would vote to cut out a little bit smaller part (maybe the lower left quarter) and to make the circles a bit larger and thicker.

OK, we can go smaller. However I would enlarge the circles as they show the actual apertures. But if you mean make them relatively larger, that having a smaller cut will solve that, too.

having a smaller cut

I don't know what you mean here.

I mean the smaller cut-out you sad. Now it's 400x400 pixels, if it's just 200x200 the fix 4pixel radius will look relatively larger for the apertures.

@cdeil @bsipocz - it looks like the apertures are systematically off the center of the sources, is there a (0.5,0.5) or (1,1) pixel offset?

Ah ... OK .. my suggestion was to make the circles so large that they visually contain the whole source, i.e. are not in the black part of the image when using the inverted color scale as we are now.

My main point would be that this makes them more visible, but I guess it also corresponds to "full-enclosure" aperture photometry, i.e. getting the total flux from each source, which is anyways a reasonable thing to do as long as we don't have a PSF model to provide a correction for the flux in smaller apertures?

@astrofrog To me it looks like some apertures are well-centered on the sources and only some are shifted ... @larrybradley maybe you can help check if the daophot results for this example and these parameters are good?

This is the version @bsipocz showed before, it seems the apertures are very slightly too low and too left. Do you see that too?

@cdeil - OK, not I see what you meant. Optimizing the aperture size for the brightest sources means on the other hand that you measure too many sky in the fainter sources, and so screw up their photometry. The brightest sources in this image are saturated anyway, so most people would throw them away. One typical approach is to use multiple apertures per sources and choose an optimal one for each in a later step.

@astrofrog - yes, I had the same feeling, but couldn't find where it went wrong.

I will compare daofind with iraf's version, but my quick reaction is that this could be an issue in matplotlib (I've seen similar plotting offsets before). I was going to look into that, but stopped fearing it would lead me down a deep rabbit hole!

@larrybradley - now that you said it, I remember we indeed had some similar issues with matplotlib.

Actually I think the issue is that matplotlib is assuming te center of the bottom left pixel is 0.5,0.5, so I think you need e.g.:

plt.imshow(image, cmap='gray_r', origin='lower', extent=[-0.5, nx - 0.5, -0.5, ny - 0.5])

Blinking the two images there's an obvious shift in the displayed image and I think there shouldn't be any.

I don't see a shift in either of these images ... both look correct to me (blinking in Chrome).
To not get sidetracked further here I created a new issue about the possible shift when plotting circles over an image: astrofrog/wcsaxes#87

Sorry, my previous comment was a bit fuzzy.
The shift I see is in the background image and only see it when the interpolation is on. It may not be even half a pixel and some of it may due to be an illusion only. I don't think there are any problems with the circle plotting.

Looking at the current example I'm not sure adding all the code to find and plot the faintest and brightest source is an improvement ... @bsipocz what do you think?

Possible alternatives are to end with the line running the aperture photometry and to simply state that the returned flux is an array. Or maybe a scatter plot of flux estimated by daophot against flux from the aperture photometry could be a very simple and quick visualisation if we want a second plot. (I didn't look what daophot does, so I'm not sure if this is meaningful).

In any case I'd like to merge this soon so that we can give Perry this plot for the presentation tomorrow ... the code will change in the next weeks anyways when aperture_photometry starts to return Table things will become simpler.

I'm not sure I like this plot. As expected we don't do well with the saturated stars:

Agreed ... let's not show the fluxes for now ... it will take more time to get reliable results.

Sorry I missed the discussion about the flux comparison, but I would not use this plot. Note, that this is daofind not DAOPHOT (which does a proper PSF fitting). daofind's flux is a bit arbitrary and depends on the kernel FWHM (I merely replicated what it does in iraf for consistency) -- I would never use to measure actual fluxes. Perhaps we could compare the aperture photometry with actual PSF-fitting photometry from photutils.psf?

@larrybradley I agree eventually we probably want to add aperture photometry to this example, but the way it's implemented at the moment it probably won't work for so few stars (see discussion in #104) and it's mostly untested, so I'd suggest leaving it out for today.

As I said somewhere above, IMO it would be OK to end the example with flux = aperture_photometry(...) for today and add the other things in the future.

I just now sent this slide to Perry for tomorrow's presentation of Astropy at Scipy where photutils will be mentioned:

Fitting the code example on a slide would have been difficult in any case. So I don't think we have to rush this PR too much ... if it's ready by tomorrow fine, but if it takes a few more days, no problem.

@cdeil Sounds good to me. (I think you meant we should eventually add PSF photometry as an example).

And thanks for putting the slide together for Scipy.

This bkg estimation is rather ugly, but for the time being it's better than nothing. Most of the code is under refactoring anyway.

@bsipocz Thanks for all the work you've put into this and the patience with my comments.

I'm not sure hard-wiring this one form of background estimation into aperture_photometry like you do here. If / how to add background estimation to photutils should first be discussed and then done in a separate pull request (with docs / tests added).

For this example I'd suggest to subtract the global median from the image at the top and add a sentence or two explaining that from here on for this example we assume the image is background-subtracted. In the future better background estimation methods might be added to photutils or ccdproc.

Before I merge this, could you please squash the commits into a single commit "Update getting started example"? (using git rebase master --interactive from this branch)

@cdeil - I don't like the hard wired way of estimating the local background either, and definitely would be good to have some astropy/ccdproc/photutils tools to do local bkg estimations.
One approach is already in the aperture_photomery docs, but it doesn't work nicely in crowed fields, like the one in this example.

@cdeil - I think this is ready to be merged. I've removed the duplicate plot, and only kept the second one.

Here's a comparison of photutils.daofind() and IRAF DAOFIND for this example, which shows excellent agreement:
http://nbviewer.ipython.org/gist/larrybradley/50666b1ba864020be8c8

Note I used photutils.utils.img_stats() to estimate the background and rms using simple sigma-clipping. I lowered the threshold here to 2.5 sigma because the thresholding is actually performed on a convolved image. The non-detected sources are likely saturated and don't have Gaussian profiles (or don't pass the sharpness/roundness criteria). Actually, it seems that many sources in this image (not just the ones with visible diffraction spikes) are nearly saturated with flat-top profiles, which may cause detection problems. In any case, photutils.daofind() and IRAF's DAOFIND are consistent.

Looks good .. merging this now.

Everyone, feel free to extend / improve this example and other parts of the docs in future pull requests.

@cdeil - I reported the bug in the shift here: matplotlib/matplotlib#3776

Actually my MWE was showing a different bug, so investigating some more.

Ok, I'm clearly confused - re-reading this I realize that there was actually no bug and that it was just previous versions that had a bug. But everything seems fine here.