code: def interp(self, energies, dtheta, scale_fn=None):
          
  
          if scale_fn is None and self.scale_fn:
              scale_fn = self.scale_fn
  
          log_energies = np.log10(energies)
  
          shape = (energies * dtheta).shape
          scale_factor = np.ones(shape)
  
          if scale_fn is not None:
              dtheta = dtheta / scale_fn(energies)
              scale_factor = 1. / scale_fn(energies)**2
  
          vals = np.exp(self._psf_fn((dtheta, log_energies)))
          return vals * scale_factor
docs: Evaluate the PSF model at an array of energies and angular
          separations.
  
          Parameters
          ----------
          energies : array_like
              Array of energies in MeV.
  
          dtheta : array_like
              Array of angular separations in degrees.
  
          scale_fn : callable        
              Function that evaluates the PSF scaling function.
              Argument is energy in MeV.