672 origVariance=
False, nptsFull=1e6, keepPlots=
True, titleFs=14):
673 """Plot diffim residuals for LOCAL and SPATIAL models.
679 for cell
in kernelCellSet.getCellList():
680 for cand
in cell.begin(
True):
682 if not (cand.getStatus() == afwMath.SpatialCellCandidate.GOOD):
685 diffim = cand.getDifferenceImage(diffimLib.KernelCandidateF.ORIG)
686 orig = cand.getScienceMaskedImage()
688 ski = afwImage.ImageD(kernel.getDimensions())
689 kernel.computeImage(ski,
False, int(cand.getXCenter()), int(cand.getYCenter()))
691 sbg =
background(int(cand.getXCenter()), int(cand.getYCenter()))
692 sdiffim = cand.getDifferenceImage(sk, sbg)
695 bbox = kernel.shrinkBBox(diffim.getBBox())
696 tdiffim = diffim.Factory(diffim, bbox)
697 torig = orig.Factory(orig, bbox)
698 tsdiffim = sdiffim.Factory(sdiffim, bbox)
701 candidateResids.append(np.ravel(tdiffim.getImage().getArray()
702 / np.sqrt(torig.getVariance().getArray())))
703 spatialResids.append(np.ravel(tsdiffim.getImage().getArray()
704 / np.sqrt(torig.getVariance().getArray())))
706 candidateResids.append(np.ravel(tdiffim.getImage().getArray()
707 / np.sqrt(tdiffim.getVariance().getArray())))
708 spatialResids.append(np.ravel(tsdiffim.getImage().getArray()
709 / np.sqrt(tsdiffim.getVariance().getArray())))
711 fullIm = diffExposure.getMaskedImage().getImage().getArray()
712 fullMask = diffExposure.getMaskedImage().getMask().getArray()
714 fullVar = exposure.getMaskedImage().getVariance().getArray()
716 fullVar = diffExposure.getMaskedImage().getVariance().getArray()
719 bitmaskBad |= afwImage.Mask.getPlaneBitMask(
'NO_DATA')
720 bitmaskBad |= afwImage.Mask.getPlaneBitMask(
'SAT')
721 idx = np.where((fullMask & bitmaskBad) == 0)
722 stride = int(len(idx[0])//nptsFull)
723 sidx = idx[0][::stride], idx[1][::stride]
724 allResids = fullIm[sidx]/np.sqrt(fullVar[sidx])
726 testFootprints = diffimTools.sourceToFootprintList(testSources, warpedTemplateExposure,
727 exposure, config, Log.getDefaultLogger())
728 for fp
in testFootprints:
729 subexp = diffExposure.Factory(diffExposure, fp[
"footprint"].getBBox())
730 subim = subexp.getMaskedImage().getImage()
732 subvar = afwImage.ExposureF(exposure, fp[
"footprint"].getBBox()).getMaskedImage().getVariance()
734 subvar = subexp.getMaskedImage().getVariance()
735 nonfitResids.append(np.ravel(subim.getArray()/np.sqrt(subvar.getArray())))
737 candidateResids = np.ravel(np.array(candidateResids))
738 spatialResids = np.ravel(np.array(spatialResids))
739 nonfitResids = np.ravel(np.array(nonfitResids))
743 from matplotlib.font_manager
import FontProperties
744 except ImportError
as e:
745 print(
"Unable to import pylab: %s" % e)
751 fig.canvas._tkcanvas._root().lift()
755 fig.suptitle(
"Diffim residuals: Normalized by sqrt(input variance)", fontsize=titleFs)
757 fig.suptitle(
"Diffim residuals: Normalized by sqrt(diffim variance)", fontsize=titleFs)
759 sp1 = pylab.subplot(221)
760 sp2 = pylab.subplot(222, sharex=sp1, sharey=sp1)
761 sp3 = pylab.subplot(223, sharex=sp1, sharey=sp1)
762 sp4 = pylab.subplot(224, sharex=sp1, sharey=sp1)
763 xs = np.arange(-5, 5.05, 0.1)
764 ys = 1./np.sqrt(2*np.pi)*np.exp(-0.5*xs**2)
766 sp1.hist(candidateResids, bins=xs, normed=
True, alpha=0.5, label=
"N(%.2f, %.2f)"
767 % (np.mean(candidateResids), np.var(candidateResids)))
768 sp1.plot(xs, ys,
"r-", lw=2, label=
"N(0,1)")
769 sp1.set_title(
"Candidates: basis fit", fontsize=titleFs - 2)
770 sp1.legend(loc=1, fancybox=
True, shadow=
True, prop=FontProperties(size=titleFs - 6))
772 sp2.hist(spatialResids, bins=xs, normed=
True, alpha=0.5, label=
"N(%.2f, %.2f)"
773 % (np.mean(spatialResids), np.var(spatialResids)))
774 sp2.plot(xs, ys,
"r-", lw=2, label=
"N(0,1)")
775 sp2.set_title(
"Candidates: spatial fit", fontsize=titleFs - 2)
776 sp2.legend(loc=1, fancybox=
True, shadow=
True, prop=FontProperties(size=titleFs - 6))
778 sp3.hist(nonfitResids, bins=xs, normed=
True, alpha=0.5, label=
"N(%.2f, %.2f)"
779 % (np.mean(nonfitResids), np.var(nonfitResids)))
780 sp3.plot(xs, ys,
"r-", lw=2, label=
"N(0,1)")
781 sp3.set_title(
"Control sample: spatial fit", fontsize=titleFs - 2)
782 sp3.legend(loc=1, fancybox=
True, shadow=
True, prop=FontProperties(size=titleFs - 6))
784 sp4.hist(allResids, bins=xs, normed=
True, alpha=0.5, label=
"N(%.2f, %.2f)"
785 % (np.mean(allResids), np.var(allResids)))
786 sp4.plot(xs, ys,
"r-", lw=2, label=
"N(0,1)")
787 sp4.set_title(
"Full image (subsampled)", fontsize=titleFs - 2)
788 sp4.legend(loc=1, fancybox=
True, shadow=
True, prop=FontProperties(size=titleFs - 6))
790 pylab.setp(sp1.get_xticklabels() + sp1.get_yticklabels(), fontsize=titleFs - 4)
791 pylab.setp(sp2.get_xticklabels() + sp2.get_yticklabels(), fontsize=titleFs - 4)
792 pylab.setp(sp3.get_xticklabels() + sp3.get_yticklabels(), fontsize=titleFs - 4)
793 pylab.setp(sp4.get_xticklabels() + sp4.get_yticklabels(), fontsize=titleFs - 4)
800 if keepPlots
and not keptPlots:
803 print(
"%s: Please close plots when done." % __name__)
808 print(
"Plots closed, exiting...")
810 atexit.register(show)
A kernel created from an Image.