docs: def create_header(coord, radius, proj='ZEA', npix=30):
      """ Create a header a new image """
      gal = coord.name == 'galactic'
      values = [
              ["NAXIS",  2,          ],
   
              ["NAXIS1", npix,       ],
              ["NAXIS2", npix,       ],
   
              ["CTYPE1", 'GLON-%s'%proj if gal else 'RA---%s'%proj ],
              ["CTYPE2", 'GLAT-%s'%proj if gal else 'DEC--%s'%proj ],
   
              ["CRPIX1", npix/2. + 0.5,       ],
              ["CRPIX2", npix/2. + 0.5,       ],
   
              ["CRVAL1", coord.l.deg if gal else coord.ra.deg,        ],
              ["CRVAL2", coord.b.deg if gal else coord.dec.deg,       ],
   
              ["CDELT1", -3.*radius/npix,       ],
              ["CDELT2",  3.*radius/npix,        ],
      ]
   
      if not gal:
          values += [
              ['RADECSYS','FK5'],
              ['EQUINOX',2000],
          ]
   
      cards = [pyfits.Card(*i) for i in values]
      header=pyfits.Header(cards=cards)
   
      return header