import poppy 
import astropy
import cloud
import matplotlib.pyplot as plt










def psf(wavelength):
    osys = poppy.OpticalSystem(oversample=10)
    osys.addPupil(function='Circle', radius=.25,pad_factor = 1.0) # pupil radius in meters
    osys.addDetector(pixelscale=0.18, fov_arcsec=23)  # image plane coordinates in arcseconds
    return osys.calcPSF(wavelength,display_intermediates=False)




wavelength=0.633e-6
manywaves=400e-9 + 1e-8*np.arange(10,41) #400 to 800  nm
x=cloud.call(psf,wavelength,_env='poppy',_profile=True)
multiwave=cloud.result(cloud.map(psf,manywaves,_env='poppy',_profile=True))

poppy.display_PSF(cloud.result(x))

#add up the returned psfs for each wavlength:
summed=np.sum([(wave[0].data) for wave in multiwave[:]],axis=0)
   
imshow(log10(summed))

hdufieldresult = astropy.io.fits.HDUList(astropy.io.fits.PrimaryHDU(np.array(summed)))
hdufieldresult.writeto('summed.fits')