Visualization utilities

[svank]

Here's a more substantive PR, involving design decisions---I expect comments and opinions! :)

I've rigged up find_stars_and_average to store the number of stars contributing to each patch, and to store those as a .counts attribute on the PatchCollection.

I added a visualization module with two functions---one to plot that counts data:

plt.figure(figsize=(7,7))
regularizepsf.visualize.visualize_patch_counts(patch_collection, ax=plt.gca())

For this plot, at first I was making an array containing the counts and imshowing that, to make a plot like this (color-by-numbers!):

--------------------
| 7 | 12 | 11 | 5  |
| 6 | 8  | 19 | 10 |
--------------------

But the patches overlap each other and stars contribute to multiple patches, so I thought that sort of plot with space-filling array cells makes it really easy to mis-interpret the plot as marking the number of stars contained within each cell. So instead I did this scatter-plot approach.

The second function plots the patches themselves. (You must also have code to do this---maybe you should substitute yours in here)

regularizepsf.visualize.visualize_patches(cpc, psf_size=psf_size, figsize=(12,12), uniform_scaling=True)

Of note, none of this has tests currently

[codecov-commenter]

Codecov Report

Patch coverage: 93.95% and project coverage change: +0.54 🎉

Comparison is base (0011641) 94.03% compared to head (02171b1) 94.57%.

📣 This organization is not using Codecov’s GitHub App Integration. We recommend you install it so Codecov can continue to function properly for your repositories. Learn more

Additional details and impacted files

@@            Coverage Diff             @@
##             main      #17      +/-   ##
==========================================
+ Coverage   94.03%   94.57%   +0.54%     
==========================================
  Files           9       11       +2     
  Lines         821      996     +175     
==========================================
+ Hits          772      942     +170     
- Misses         49       54       +5     

Impacted Files	Coverage Δ
setup.py	0.00% <ø> (ø)
regularizepsf/visualize.py	91.45% <91.45%> (ø)
tests/test_visualize.py	98.00% <98.00%> (ø)
regularizepsf/__init__.py	100.00% <100.00%> (ø)
regularizepsf/fitter.py	90.00% <100.00%> (+3.30%)	⬆️

Help us with your feedback. Take ten seconds to tell us how you rate us. Have a feature suggestion? Share it here.

☔ View full report in Codecov by Sentry.
📢 Do you have feedback about the report comment? Let us know in this issue.

[jmbhughes]

I haven't looked at your code yet, but I will by the end of Monday.

As you suggested, I might merge some of my patch plotting code into yours. I might also add a function to plot both a before and after since that's a comparison that's nice to have.

Annoyingly, this plot has a transparent background that doesn't render well in GitHub's dark mode.

As for tests, I'm concerned how this pulls down test coverage. Writing good unit tests for Matplotlib code can be tricky. A quick Google search gives many suggestions on how to unit test figure generation code. Even just testing that it runs would be better than nothing.

[jmbhughes]

Also, I think we need to add matplotlib to the setup.py install_requires list. Probably should validate that contains everything in requirements.txt too. This makes me wonder if the CI should be using that list instead of requirements.txt as well.

[svank]

Your patch plot looks better than mine.

Annoyingly, this plot has a transparent background that doesn't render well in GitHub's dark mode.

Jupyter? If you drop c.InlineBackend.rc = {'figure.facecolor': 'white'} in ~/.ipython/profile_default/ipython_kernel_config.py, matplotlib will fill plots in with a solid white background (and I wish this were the default!)

I have some matplotlib-based tests in https://github.com/svank/wispr_analysis, but they feel really flakey---especially anything that includes text, since apparently font rendering is incredibly platform-specific

[jmbhughes]

We should make a page in the Jupyter-book docs website for this module.

[jmbhughes]

Why would you not always want uniform_scaling to be True?

An alternative way to do this would be to allow specifying a vmin and vmax. If they're set to None then it lets that value vary by different patches. If they're specified, then all the patches use that vmin and vmax.

[jmbhughes]

Maybe this should be specified in the imshow_args just to avoid redundant variables. An example could show how to do it.

[svank]

I think the data range should generally be [0, 1] for every patch, since the star cutouts are each normalized, but I had been seeing some of my patches span that range and others only go up to ~0.6 or so, so I thought there might be use cases for either showing how each patch compares to the others, or using the full dynamic range for each individual patch plot.

But I think that variation was because of the way the identified stellar coordinates were being truncated rather than rounded, so the peak-value pixels weren't always aligned. After changing the code to round the stellar coordinates, I think I'm seeing nice [0, 1] ranges for each patch, so this may be unnecessary now. I think we can yank out the uniform_scaling arg, have the default imshow_args set vmin=0, vmax=1, and the user can override those to None if they really want to.

[jmbhughes]

I'm glad you allow imshow_args because specifying a custom colormap is important. For the paper, we used a gamma-corrected colormap. We might want to supply that as an option for easy use. This is how I define the paper's colormap:

import matplotlib.colors as clr

a = np.linspace(0, 1, 1000)
r = np.sqrt(a)
g = a
b = np.square(a)
colors = np.stack([r, g, b], axis=-1)
custom = clr.ListedColormap(colors)

We also specified a norm option for the imshow command:

  axs[i, j].imshow(array_corrector[x, y][117:-117, 117:-117], 
                   cmap=custom, 
                   norm=colors.PowerNorm(gamma=1/2.2, vmin=vmin, vmax=vmax))

It would be good to call out both in the documentation website and in the docstring here the norm and cmap parts of this as good practice.

[svank]

Ah, I was looking through the default matplotlib colormaps and couldn't figure out which one you were using---this is why! And the gamma correction is probably a good default

[jmbhughes]

Can we allow a user to specify a custom trim? Maybe they specify a large PSF size for some reason but want the plot just to show the center pixels in some box.

[svank]

Good call!

[jmbhughes]

I think we have different approaches in our plotting here. Since each pixel is sampled 4 times, there are essentially 4 grids of patches. Each grid completely covers the image without overlapping other patches in that grid. You can see how the grids are constructed in this part of the code. In my plots, I took the first quarter of patches and plotted them since there's no overlap. It looks like you're using all the patches. We should discuss which approach is best.

[svank]

I think it should be possible both to plot a fraction of the patches, and to plot all the patches. The former may be good for cleaner illustrations, while the latter is good for checking that all the patches look reasonable.

One thing I was unsure about in my code was plotting all the patches (from all four grids) together inside one large set of axes labeled with pixel numbers. The axis labels help connect a patch to an image location, but they misrepresent the range of each patch, since I was cramming them all into one big plot despite their sample overlapping regions. I like that your approach removes all the axis labels and adds small borders between each patch, to make their locations more representative than pixel-perfect.

[jmbhughes]

Here's the code I used to visualize a before and after set of patches:

import matplotlib.pyplot as plt
from matplotlib.gridspec import GridSpec
from mpl_toolkits.axes_grid1.inset_locator import inset_axes
import itertools
import matplotlib.pylab as pylab
import matplotlib.colors as colors
params = {'axes.titlesize': 15}
pylab.rcParams.update(params)

vmin = 0.02
vmax = 1.0

fig = plt.figure(constrained_layout=True, figsize=(13, 6))

gs = GridSpec(4, 9, figure=fig, 
              width_ratios=([1]*4) + [0.8] + ([1]*4),
              height_ratios=[1]*4)

uncorrected_axs = dict()
for i, j in itertools.product(range(4), range(4)):
    uncorrected_axs[i, j] = fig.add_subplot(gs[i, j])
    
for i in range(16):
    row, col = i % 4, i // 4
    x = coords[np.array([1, 3, 5, 7])][row]
    y = coords[np.array([1, 3, 5, 7])][col]
    im = uncorrected_axs[3-row, col].imshow(array_corrector[x, y][117:-117, 117:-117], 
                                  origin='lower', cmap=custom, 
                                            norm=colors.PowerNorm(gamma=1/2.2, vmin=vmin, vmax=vmax))
    uncorrected_axs[3-row, col].set_axis_off()
    uncorrected_axs[3-row, col].set_aspect(1)
    
cax = fig.add_axes([0.479, 0.065, 0.012, 0.85])
#fig.colorbar(im, cax=cax, orientation='vertical', label='log$_{10}$(brightness)')
fig.colorbar(im, cax=cax, orientation='vertical', label='brightness')

corrected_axs = dict()
for i, j in itertools.product(range(4), range(4)):
    corrected_axs[i, j] = fig.add_subplot(gs[i, j+5])
for i in range(16):
    row, col = i // 4, i % 4
    x = coords[np.array([1, 3, 5, 7])][row]
    y = coords[np.array([1, 3, 5, 7])][col]
    im = corrected_axs[3-row, col].imshow(output[CoordinateIdentifier(None, x, y)][117:-117, 117:-117], 
                                    origin='lower', cmap=custom, 
                                          norm=colors.PowerNorm(gamma=1/2.2, vmin=vmin, vmax=vmax))
    corrected_axs[3-row, col].set_axis_off()
    corrected_axs[3-row, col].set_aspect(1)

    
plt.figtext(0.23, 0.95, 'Uncorrected', ha='center', va='center', fontsize=15)
plt.figtext(0.23+0.54, 0.95, 'Corrected', ha='center', va='center', fontsize=15)
uncorrected_axs[0, 2].set_title("\n")
plt.show()
fig.savefig("/Users/jhughes/Desktop/punch.png", dpi=300)

It would obviously need generalizing if we were to add it.

[svank]

It would obviously need generalizing if we were to add it.

fig.savefig("/Users/jhughes/Desktop/punch.png", dpi=300)

I think we can just require the user to log in as jhughes, no need to generalize

[svank]

But I do prefer your patch plots, so I think this code should be the base of what we end up on

[jmbhughes]

I think the type hint of counts should be Optional since we're setting it to None.

[jmbhughes]

As for testing, it's not absolutely necessary. I think we can just say in any software reviews that testing plotting is finicky. 88% isn't bad. (Maybe we could boost it a little by checking if there's other code not covered.)

It looks like you can specify to not show text when running a Matplotlib test. Maybe that would help?

[svank]

I've updated to roll in your patch-plotting code.

I have a before-and-after mode (in this demo the "corrected" side isn't actually corrected, because I don't have that ready yet):
regularizepsf.visualize.visualize_PSFs(array_corrector, cpc, False, region_size=psf_size)

and a plain "what's my PSF" mode:
regularizepsf.visualize.visualize_PSFs(array_corrector, None, False, region_size=psf_size)

Both modes can also plot all the patches:
regularizepsf.visualize.visualize_PSFs(array_corrector, None, True, region_size=psf_size)

(Docs and tests are still TODO)

[svank]

I've changed visualize_patch_counts to imshow an array of count values, rather than make a scatter plot. It makes for a cleaner plot, I think, and doesn't result in overlapping scatterplot dots when there's a large number of patches. My previous reason for not doing it this way was that the axis labels didn't really make sense as pixel numbers, but I've changed it to label the axes as just "patch number" instead.

I've also added first-draft tests for these two functions. I used pytest-mpl, which has the advantage over the matplotlib.testing tools that it doesn't save out the generated plots when the tests pass. It doesn't run the image comparison tests at all unless you run the tests as pytest --mpl, and there's another flag for re-generating the reference images, so I added a stub "Development" page to the docs describing this. I think having the image tests not run by default is so that you can restrict them to a consistent test environment with consistent font rendering. That can be worked around by removing all text from the plots---I was disappointed that the remove_text option only removes some text, so I went heavy on removing all text from the image tests.

[svank]

Just continuing to add more for you to look at when you can, I added a tool to compute and plot the transfer kernel corresponding to each patch. Looks like this:
regularizepsf.visualize.visualize_transfer_kernels(array_corrector, alpha, epsilon, region_size=psf_size)

To make this work easily, I went into the Cython code and pulled the regularized inverse into its own function.

Something that occurred to me is that for the image-comparison tests I added, the reference images are generated with the changes of #19 and #28 in effect, so the tests won't pass as-is unless this PR is applied on top of those.

[svank]

Another update I'm proud of: labeling the axes of these plots with the pixel bounds of each patch. I found this really helpful for matching up these plots to the actual data. It makes the axes a bit more verbose, though, so I have this feature default to off---for large images or small patch sizes, there could be ugly overlapping text.

I tried two styles, and should probably choose just one for everything. The second is a bit cleaner but takes up more space.

[jmbhughes]

This looks good. I think I like the second one better too.

[jmbhughes]

@svank I went ahead and merged the other PRs.

[svank]

I've added docs---I think this PR is in a mergeable state now

[jmbhughes]

Awesome! Thanks for adding all of this!