aboutsummaryrefslogtreecommitdiffstats
path: root/bin/hpx2healpix
blob: ce83553d616d8c712860013e61c24904fc4c1650 (plain)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
#!/usr/bin/env python3
#
# Copyright (c) 2016 Weitian LI <liweitianux@live.com>
# MIT license

"""
Recover the sky map in HPX projection format back into HEALPix table format,
i.e., the reverse of `healpix2hpx.py`.
"""

import os
import sys
import argparse
import logging

import numpy as np
from astropy.io import fits

import fg21sim
from fg21sim.configs import configs
from fg21sim.utils import hpx2healpix, setup_logging


# Reference:
# http://docs.astropy.org/en/stable/io/fits/usage/table.html#column-creation
FITS_COLUMN_FORMATS = {
    np.dtype("int16"): "I",
    np.dtype("int32"): "J",
    np.dtype("int64"): "K",
    np.dtype("float32"): "E",
    np.dtype("float64"): "D",
}


def main():
    parser = argparse.ArgumentParser(
        description="Recover the image in HPX projection to HEALPix data")
    parser.add_argument("infile", help="input FITS image in HPX projection")
    parser.add_argument("outfile", help="output HEALPix data file")
    parser.add_argument("-C", "--clobber", action="store_true",
                        help="overwrite the existing output file")
    parser.add_argument("-F", "--float", action="store_true",
                        help="use float (single precision) instead of double")
    parser.add_argument("-l", "--log", dest="loglevel", default=None,
                        choices=["DEBUG", "INFO", "WARNING",
                                 "ERROR", "CRITICAL"],
                        help="set the log level")
    parser.add_argument("-L", "--logfile", default=None,
                        help="filename where to save the log messages")
    parser.add_argument("-Q", "--quiet", action="store_true",
                        help="be quiet so do not log messages to screen")
    args = parser.parse_args()

    if args.quiet:
        log_stream = ""
    else:
        log_stream = None

    tool = os.path.basename(sys.argv[0])
    pkgname = fg21sim.__pkgname__

    setup_logging(dict_config=configs.logging,
                  level=args.loglevel,
                  stream=log_stream,
                  logfile=args.logfile)
    logger = logging.getLogger(tool)
    logger.info("COMMAND: {0}".format(" ".join(sys.argv)))

    history = [
        "TOOL: {0}".format(tool),
        "PARAM: {0}".format(" ".join(sys.argv[1:])),
    ]
    comments = [
        'Tool "{0}" is part of the "{1}" package'.format(tool, pkgname),
        'distributed under {0} license.'.format(fg21sim.__license__),
        'See also {0}'.format(fg21sim.__url__)
    ]

    hp_data, hp_header = hpx2healpix(args.infile,
                                     append_history=history,
                                     append_comment=comments)
    if args.float:
        logger.info("HEALPix data: use single-precision float numbers")
        hp_data = hp_data.astype(np.float32)
    hdu = fits.BinTableHDU.from_columns([
        fits.Column(name="I", array=hp_data,
                    format=FITS_COLUMN_FORMATS.get(hp_data.dtype))
    ], header=hp_header)
    hdu.writeto(args.outfile, clobber=args.clobber, checksum=True)
    logger.info("HEALPix data write to FITS file: %s" % args.outfile)


if __name__ == "__main__":
    main()