GammaLib

# Create allsky images of solid angle for various projections

import gammalib

#projection = 'AIT'

#projection = 'CAR'

#projection = 'TAN'

#projection = 'MER'

#projection = 'AZP'

projection = 'STG'

x, y = 0, 0

#dx, dy = 1, 1

#nx, ny = 360, 180

dx, dy = 0.5, 0.5

nx, ny = 720, 380

nmaps = 1

image = gammalib.GSkymap(projection, "CEL", x, y,

                         -dx, dy, nx, ny, nmaps)

# Fill the sky map with the model image

total = 0.0

for pix in range(image.npix()):

    try:

        solidangle = image.solidangle(pix)

        total      = total + solidangle

        image[pix] = gammalib.rad2deg/dx * gammalib.rad2deg/dy * solidangle

    except RuntimeError:

        pass

print('Total solid angle: %f sr' % total)

# Save the image to a FITS file

filename = 'test_{0}.fits'.format(projection)

print('Writing {0}'.format(filename))

image.save(filename, True)