import matplotlib.pyplot as plt
from astropy.visualization import simple_norm
from lacosmic.utils import make_cosmic_rays, make_gaussian_sources

shape = (512, 512)
data, error = make_gaussian_sources(shape, seed=0)
norm = simple_norm(data, 'sqrt', percent=99.5)
cr_img = make_cosmic_rays(shape, n_cosmics=200, seed=0)
data2 = data + cr_img

fig, ax = plt.subplots()
axim = ax.imshow(data2, norm=norm, origin='lower')
plt.tight_layout()
plt.show()