def containment_radius(self, energy, theta, fraction=0.68):
    'Compute containment for all energy and theta values'
    energies = u.Quantity(energy).flatten()
    thetas = Angle(theta).flatten()
    radius = np.empty((theta.size, energy.size))
    for (idx, energy) in enumerate(energies):
        for (jdx, theta) in enumerate(thetas):
            try:
                psf = self.psf_at_energy_and_theta(energy, theta)
                radius[(jdx, idx)] = psf.containment_radius(fraction)
            except ValueError:
                log.debug(f'Computing containment failed for energy = {energy:.2f} and theta={theta:.2f}')
                log.debug(f'Sigmas: {psf.sigmas} Norms: {psf.norms}')
                radius[(jdx, idx)] = np.nan
    return Angle(radius, 'deg')
Compute containment for all energy and theta values
def containment_radius(self, energy, theta, fraction=0.68):
    
    energies = u.Quantity(energy).flatten()
    thetas = Angle(theta).flatten()
    radius = np.empty((theta.size, energy.size))
    for (idx, energy) in enumerate(energies):
        for (jdx, theta) in enumerate(thetas):
            try:
                psf = self.psf_at_energy_and_theta(energy, theta)
                radius[(jdx, idx)] = psf.containment_radius(fraction)
            except ValueError:
                log.debug(f'Computing containment failed for energy = {energy:.2f} and theta={theta:.2f}')
                log.debug(f'Sigmas: {psf.sigmas} Norms: {psf.norms}')
                radius[(jdx, idx)] = np.nan
    return Angle(radius, 'deg')