Merge pull request #6 from Andreyousan/code_trimming

Code trimming
This commit is contained in:
Андрей Мухин 2022-12-15 15:35:37 +03:00 committed by GitHub
commit a4698b3bee
No known key found for this signature in database
GPG Key ID: 4AEE18F83AFDEB23
21 changed files with 613 additions and 79 deletions

View File

@ -1,2 +1,3 @@
include nuwavsource/pixpos/nuApixpos20100101v007.fits
include nuwavsource/pixpos/nuBpixpos20100101v007.fits
include nuwavsource/pixpos/nuBpixpos20100101v007.fits
include nuwavsource/badpix_headers/*

View File

@ -0,0 +1,46 @@
XTENSION= 'BINTABLE' / binary table extension
BITPIX = 8 / 8-bit bytes
NAXIS = 2 / 2-dimensional binary table
NAXIS1 = 12 / width of table in bytes
NAXIS2 = 64 / number of rows in table
PCOUNT = 0 / size of special data area
GCOUNT = 1 / one data group (required keyword)
TFIELDS = 4 / number of fields in each row
TTYPE1 = 'TIME ' / Start Time of Bad Pixel Interval
TFORM1 = '1D ' / data format of field: 8-byte DOUBLE
TUNIT1 = 's ' / physical unit of field
TTYPE2 = 'RAWX ' / X-position in Raw Detector Coordinates
TFORM2 = '1B ' / data format of field: BYTE
TUNIT2 = 'pixel ' / physical unit of field
TTYPE3 = 'RAWY ' / Y-position in Raw Detector Coordinates
TFORM3 = '1B ' / data format of field: BYTE
TUNIT3 = 'pixel ' / physical unit of field
TTYPE4 = 'BADFLAG ' / Bad pixel flag
TFORM4 = '16X ' / data format of field: BIT
EXTNAME = 'BADPIX ' / Name of the binary table extension
EXTVER = 1 / There shall be one instance of this extension f
DETNAM = 'DET0 ' / CZT Detector ID (0,1,2 or 3)
TELESCOP= 'NuSTAR ' / Telescope (mission) name
INSTRUME= 'FPMA ' / Instrument name (FPMA or FPMB)
ORIGIN = 'Caltech ' / Source of FITS file
CREATOR = 'FTOOLS 6.10 ' / Creator
VERSION = 1 / Extension version number
FILENAME= 'nuAuserbadpix20100101v001.fits' / File name
CONTENT = 'NuSTAR User Bad Pixel Table' / File content
CCLS0001= 'BCF ' / Daset is a Basic Calibration File
CDTP0001= 'DATA ' / Calibration file contains data
CCNM0001= 'BADPIX ' / Type of calibration data
CVSD0001= '2010-01-01' / Date when this file should first be used
CVST0001= '00:00:00' / Time of day when this file should first be used
CDES0001= 'NuSTAR User Bad Pixel Table' / Description
COMMENT
COMMENT This extension provides, for Detector #0 of the FPMA, the list of
COMMENT user-defined bad pixels.
COMMENT The BADFLAG column contains the bad pixel flags. These are a
COMMENT combination of the following bit settings:
COMMENT
COMMENT b0000000000000001 Bad pixel from on-ground CALDB Bad Pixel File
COMMENT b0000000000000010 Disabled pixel from on-board software
COMMENT b0000000000000100 Bad pixels in the file provided by the user
COMMENT
END

View File

@ -0,0 +1,46 @@
XTENSION= 'BINTABLE' / binary table extension
BITPIX = 8 / 8-bit bytes
NAXIS = 2 / 2-dimensional binary table
NAXIS1 = 12 / width of table in bytes
NAXIS2 = 64 / number of rows in table
PCOUNT = 0 / size of special data area
GCOUNT = 1 / one data group (required keyword)
TFIELDS = 4 / number of fields in each row
TTYPE1 = 'TIME ' / Start Time of Bad Pixel Interval
TFORM1 = '1D ' / data format of field: 8-byte DOUBLE
TUNIT1 = 's ' / physical unit of field
TTYPE2 = 'RAWX ' / X-position in Raw Detector Coordinates
TFORM2 = '1B ' / data format of field: BYTE
TUNIT2 = 'pixel ' / physical unit of field
TTYPE3 = 'RAWY ' / Y-position in Raw Detector Coordinates
TFORM3 = '1B ' / data format of field: BYTE
TUNIT3 = 'pixel ' / physical unit of field
TTYPE4 = 'BADFLAG ' / Bad pixel flag
TFORM4 = '16X ' / data format of field: BIT
EXTNAME = 'BADPIX ' / Name of the binary table extension
EXTVER = 2 / There shall be one instance of this extension f
DETNAM = 'DET1 ' / CZT Detector ID (0,1,2 or 3)
TELESCOP= 'NuSTAR ' / Telescope (mission) name
INSTRUME= 'FPMA ' / Instrument name (FPMA or FPMB)
ORIGIN = 'Caltech ' / Source of FITS file
CREATOR = 'FTOOLS 6.10 ' / Creator
VERSION = 1 / Extension version number
FILENAME= 'nuAuserbadpix20100101v001.fits' / File name
CONTENT = 'NuSTAR User Bad Pixel Table' / File content
CCLS0001= 'BCF ' / Daset is a Basic Calibration File
CDTP0001= 'DATA ' / Calibration file contains data
CCNM0001= 'BADPIX ' / Type of calibration data
CVSD0001= '2010-01-01' / Date when this file should first be used
CVST0001= '00:00:00' / Time of day when this file should first be used
CDES0001= 'NuSTAR User Bad Pixel Table' / Description
COMMENT
COMMENT This extension provides, for Detector #1 of the FPMA, the list of
COMMENT user-defined bad pixels.
COMMENT The BADFLAG column contains the bad pixel flags. These are a
COMMENT combination of the following bit settings:
COMMENT
COMMENT b0000000000000001 Bad pixel from on-ground CALDB Bad Pixel File
COMMENT b0000000000000010 Disabled pixel from on-board software
COMMENT b0000000000000100 Bad pixels in the file provided by the user
COMMENT
END

View File

@ -0,0 +1,46 @@
XTENSION= 'BINTABLE' / binary table extension
BITPIX = 8 / 8-bit bytes
NAXIS = 2 / 2-dimensional binary table
NAXIS1 = 12 / width of table in bytes
NAXIS2 = 64 / number of rows in table
PCOUNT = 0 / size of special data area
GCOUNT = 1 / one data group (required keyword)
TFIELDS = 4 / number of fields in each row
TTYPE1 = 'TIME ' / Start Time of Bad Pixel Interval
TFORM1 = '1D ' / data format of field: 8-byte DOUBLE
TUNIT1 = 's ' / physical unit of field
TTYPE2 = 'RAWX ' / X-position in Raw Detector Coordinates
TFORM2 = '1B ' / data format of field: BYTE
TUNIT2 = 'pixel ' / physical unit of field
TTYPE3 = 'RAWY ' / Y-position in Raw Detector Coordinates
TFORM3 = '1B ' / data format of field: BYTE
TUNIT3 = 'pixel ' / physical unit of field
TTYPE4 = 'BADFLAG ' / Bad pixel flag
TFORM4 = '16X ' / data format of field: BIT
EXTNAME = 'BADPIX ' / Name of the binary table extension
EXTVER = 3 / There shall be one instance of this extension f
DETNAM = 'DET2 ' / CZT Detector ID (0,1,2 or 3)
TELESCOP= 'NuSTAR ' / Telescope (mission) name
INSTRUME= 'FPMA ' / Instrument name (FPMA or FPMB)
ORIGIN = 'Caltech ' / Source of FITS file
CREATOR = 'FTOOLS 6.10 ' / Creator
VERSION = 1 / Extension version number
FILENAME= 'nuAuserbadpix20100101v001.fits' / File name
CONTENT = 'NuSTAR User Bad Pixel Table' / File content
CCLS0001= 'BCF ' / Daset is a Basic Calibration File
CDTP0001= 'DATA ' / Calibration file contains data
CCNM0001= 'BADPIX ' / Type of calibration data
CVSD0001= '2010-01-01' / Date when this file should first be used
CVST0001= '00:00:00' / Time of day when this file should first be used
CDES0001= 'NuSTAR User Bad Pixel Table' / Description
COMMENT
COMMENT This extension provides, for Detector #2 of the FPMA, the list of
COMMENT user-defined bad pixels.
COMMENT The BADFLAG column contains the bad pixel flags. These are a
COMMENT combination of the following bit settings:
COMMENT
COMMENT b0000000000000001 Bad pixel from on-ground CALDB Bad Pixel File
COMMENT b0000000000000010 Disabled pixel from on-board software
COMMENT b0000000000000100 Bad pixels in the file provided by the user
COMMENT
END

View File

@ -0,0 +1,46 @@
XTENSION= 'BINTABLE' / binary table extension
BITPIX = 8 / 8-bit bytes
NAXIS = 2 / 2-dimensional binary table
NAXIS1 = 12 / width of table in bytes
NAXIS2 = 64 / number of rows in table
PCOUNT = 0 / size of special data area
GCOUNT = 1 / one data group (required keyword)
TFIELDS = 4 / number of fields in each row
TTYPE1 = 'TIME ' / Start Time of Bad Pixel Interval
TFORM1 = '1D ' / data format of field: 8-byte DOUBLE
TUNIT1 = 's ' / physical unit of field
TTYPE2 = 'RAWX ' / X-position in Raw Detector Coordinates
TFORM2 = '1B ' / data format of field: BYTE
TUNIT2 = 'pixel ' / physical unit of field
TTYPE3 = 'RAWY ' / Y-position in Raw Detector Coordinates
TFORM3 = '1B ' / data format of field: BYTE
TUNIT3 = 'pixel ' / physical unit of field
TTYPE4 = 'BADFLAG ' / Bad pixel flag
TFORM4 = '16X ' / data format of field: BIT
EXTNAME = 'BADPIX ' / Name of the binary table extension
EXTVER = 4 / There shall be one instance of this extension f
DETNAM = 'DET3 ' / CZT Detector ID (0,1,2 or 3)
TELESCOP= 'NuSTAR ' / Telescope (mission) name
INSTRUME= 'FPMA ' / Instrument name (FPMA or FPMB)
ORIGIN = 'Caltech ' / Source of FITS file
CREATOR = 'FTOOLS 6.10 ' / Creator
VERSION = 1 / Extension version number
FILENAME= 'nuAuserbadpix20100101v001.fits' / File name
CONTENT = 'NuSTAR User Bad Pixel Table' / File content
CCLS0001= 'BCF ' / Daset is a Basic Calibration File
CDTP0001= 'DATA ' / Calibration file contains data
CCNM0001= 'BADPIX ' / Type of calibration data
CVSD0001= '2010-01-01' / Date when this file should first be used
CVST0001= '00:00:00' / Time of day when this file should first be used
CDES0001= 'NuSTAR User Bad Pixel Table' / Description
COMMENT
COMMENT This extension provides, for Detector #3 of the FPMA, the list of
COMMENT user-defined bad pixels.
COMMENT The BADFLAG column contains the bad pixel flags. These are a
COMMENT combination of the following bit settings:
COMMENT
COMMENT b0000000000000001 Bad pixel from on-ground CALDB Bad Pixel File
COMMENT b0000000000000010 Disabled pixel from on-board software
COMMENT b0000000000000100 Bad pixels in the file provided by the user
COMMENT
END

View File

@ -0,0 +1,11 @@
SIMPLE = T / file does conform to FITS standard
BITPIX = 16 / number of bits per data pixel
NAXIS = 0 / number of data axes
EXTEND = T / FITS dataset may contain extensions
COMMENT FITS (Flexible Image Transport System) format is defined in 'Astronomy
COMMENT and Astrophysics', volume 376, page 359; bibcode: 2001A&A...376..359H
TELESCOP= 'NuSTAR ' / Telescope (mission) name
INSTRUME= 'FPMA ' / Instrument name (FPMA or FPMB)
ORIGIN = 'Caltech ' / Source of FITS file
CREATOR = 'FTOOLS 6.10 ' / Creator
END

View File

@ -0,0 +1,46 @@
XTENSION= 'BINTABLE' / binary table extension
BITPIX = 8 / 8-bit bytes
NAXIS = 2 / 2-dimensional binary table
NAXIS1 = 12 / width of table in bytes
NAXIS2 = 64 / number of rows in table
PCOUNT = 0 / size of special data area
GCOUNT = 1 / one data group (required keyword)
TFIELDS = 4 / number of fields in each row
TTYPE1 = 'TIME ' / Start Time of Bad Pixel Interval
TFORM1 = '1D ' / data format of field: 8-byte DOUBLE
TUNIT1 = 's ' / physical unit of field
TTYPE2 = 'RAWX ' / X-position in Raw Detector Coordinates
TFORM2 = '1B ' / data format of field: BYTE
TUNIT2 = 'pixel ' / physical unit of field
TTYPE3 = 'RAWY ' / Y-position in Raw Detector Coordinates
TFORM3 = '1B ' / data format of field: BYTE
TUNIT3 = 'pixel ' / physical unit of field
TTYPE4 = 'BADFLAG ' / Bad pixel flag
TFORM4 = '16X ' / data format of field: BIT
EXTNAME = 'BADPIX ' / Name of the binary table extension
EXTVER = 1 / There shall be one instance of this extension f
DETNAM = 'DET0 ' / CZT Detector ID (0,1,2 or 3)
TELESCOP= 'NuSTAR ' / Telescope (mission) name
INSTRUME= 'FPMB ' / Instrument name (FPMA or FPMB)
ORIGIN = 'Caltech ' / Source of FITS file
CREATOR = 'FTOOLS 6.10 ' / Creator
VERSION = 1 / Extension version number
FILENAME= 'nuBuserbadpix20100101v001.fits' / File name
CONTENT = 'NuSTAR User Bad Pixel Table' / File content
CCLS0001= 'BCF ' / Daset is a Basic Calibration File
CDTP0001= 'DATA ' / Calibration file contains data
CCNM0001= 'BADPIX ' / Type of calibration data
CVSD0001= '2010-01-01' / Date when this file should first be used
CVST0001= '00:00:00' / Time of day when this file should first be used
CDES0001= 'NuSTAR User Bad Pixel Table' / Description
COMMENT
COMMENT This extension provides, for Detector #0 of the FPMB, the list of
COMMENT user-defined bad pixels.
COMMENT The BADFLAG column contains the bad pixel flags. These are a
COMMENT combination of the following bit settings:
COMMENT
COMMENT b0000000000000001 Bad pixel from on-ground CALDB Bad Pixel File
COMMENT b0000000000000010 Disabled pixel from on-board software
COMMENT b0000000000000100 Bad pixels in the file provided by the user
COMMENT
END

View File

@ -0,0 +1,46 @@
XTENSION= 'BINTABLE' / binary table extension
BITPIX = 8 / 8-bit bytes
NAXIS = 2 / 2-dimensional binary table
NAXIS1 = 12 / width of table in bytes
NAXIS2 = 64 / number of rows in table
PCOUNT = 0 / size of special data area
GCOUNT = 1 / one data group (required keyword)
TFIELDS = 4 / number of fields in each row
TTYPE1 = 'TIME ' / Start Time of Bad Pixel Interval
TFORM1 = '1D ' / data format of field: 8-byte DOUBLE
TUNIT1 = 's ' / physical unit of field
TTYPE2 = 'RAWX ' / X-position in Raw Detector Coordinates
TFORM2 = '1B ' / data format of field: BYTE
TUNIT2 = 'pixel ' / physical unit of field
TTYPE3 = 'RAWY ' / Y-position in Raw Detector Coordinates
TFORM3 = '1B ' / data format of field: BYTE
TUNIT3 = 'pixel ' / physical unit of field
TTYPE4 = 'BADFLAG ' / Bad pixel flag
TFORM4 = '16X ' / data format of field: BIT
EXTNAME = 'BADPIX ' / Name of the binary table extension
EXTVER = 2 / There shall be one instance of this extension f
DETNAM = 'DET1 ' / CZT Detector ID (0,1,2 or 3)
TELESCOP= 'NuSTAR ' / Telescope (mission) name
INSTRUME= 'FPMB ' / Instrument name (FPMA or FPMB)
ORIGIN = 'Caltech ' / Source of FITS file
CREATOR = 'FTOOLS 6.10 ' / Creator
VERSION = 1 / Extension version number
FILENAME= 'nuBuserbadpix20100101v001.fits' / File name
CONTENT = 'NuSTAR User Bad Pixel Table' / File content
CCLS0001= 'BCF ' / Daset is a Basic Calibration File
CDTP0001= 'DATA ' / Calibration file contains data
CCNM0001= 'BADPIX ' / Type of calibration data
CVSD0001= '2010-01-01' / Date when this file should first be used
CVST0001= '00:00:00' / Time of day when this file should first be used
CDES0001= 'NuSTAR User Bad Pixel Table' / Description
COMMENT
COMMENT This extension provides, for Detector #1 of the FPMA, the list of
COMMENT user-defined bad pixels.
COMMENT The BADFLAG column contains the bad pixel flags. These are a
COMMENT combination of the following bit settings:
COMMENT
COMMENT b0000000000000001 Bad pixel from on-ground CALDB Bad Pixel File
COMMENT b0000000000000010 Disabled pixel from on-board software
COMMENT b0000000000000100 Bad pixels in the file provided by the user
COMMENT
END

View File

@ -0,0 +1,46 @@
XTENSION= 'BINTABLE' / binary table extension
BITPIX = 8 / 8-bit bytes
NAXIS = 2 / 2-dimensional binary table
NAXIS1 = 12 / width of table in bytes
NAXIS2 = 64 / number of rows in table
PCOUNT = 0 / size of special data area
GCOUNT = 1 / one data group (required keyword)
TFIELDS = 4 / number of fields in each row
TTYPE1 = 'TIME ' / Start Time of Bad Pixel Interval
TFORM1 = '1D ' / data format of field: 8-byte DOUBLE
TUNIT1 = 's ' / physical unit of field
TTYPE2 = 'RAWX ' / X-position in Raw Detector Coordinates
TFORM2 = '1B ' / data format of field: BYTE
TUNIT2 = 'pixel ' / physical unit of field
TTYPE3 = 'RAWY ' / Y-position in Raw Detector Coordinates
TFORM3 = '1B ' / data format of field: BYTE
TUNIT3 = 'pixel ' / physical unit of field
TTYPE4 = 'BADFLAG ' / Bad pixel flag
TFORM4 = '16X ' / data format of field: BIT
EXTNAME = 'BADPIX ' / Name of the binary table extension
EXTVER = 3 / There shall be one instance of this extension f
DETNAM = 'DET2 ' / CZT Detector ID (0,1,2 or 3)
TELESCOP= 'NuSTAR ' / Telescope (mission) name
INSTRUME= 'FPMB ' / Instrument name (FPMA or FPMB)
ORIGIN = 'Caltech ' / Source of FITS file
CREATOR = 'FTOOLS 6.10 ' / Creator
VERSION = 1 / Extension version number
FILENAME= 'nuBuserbadpix20100101v001.fits' / File name
CONTENT = 'NuSTAR User Bad Pixel Table' / File content
CCLS0001= 'BCF ' / Daset is a Basic Calibration File
CDTP0001= 'DATA ' / Calibration file contains data
CCNM0001= 'BADPIX ' / Type of calibration data
CVSD0001= '2010-01-01' / Date when this file should first be used
CVST0001= '00:00:00' / Time of day when this file should first be used
CDES0001= 'NuSTAR User Bad Pixel Table' / Description
COMMENT
COMMENT This extension provides, for Detector #2 of the FPMA, the list of
COMMENT user-defined bad pixels.
COMMENT The BADFLAG column contains the bad pixel flags. These are a
COMMENT combination of the following bit settings:
COMMENT
COMMENT b0000000000000001 Bad pixel from on-ground CALDB Bad Pixel File
COMMENT b0000000000000010 Disabled pixel from on-board software
COMMENT b0000000000000100 Bad pixels in the file provided by the user
COMMENT
END

View File

@ -0,0 +1,46 @@
XTENSION= 'BINTABLE' / binary table extension
BITPIX = 8 / 8-bit bytes
NAXIS = 2 / 2-dimensional binary table
NAXIS1 = 12 / width of table in bytes
NAXIS2 = 64 / number of rows in table
PCOUNT = 0 / size of special data area
GCOUNT = 1 / one data group (required keyword)
TFIELDS = 4 / number of fields in each row
TTYPE1 = 'TIME ' / Start Time of Bad Pixel Interval
TFORM1 = '1D ' / data format of field: 8-byte DOUBLE
TUNIT1 = 's ' / physical unit of field
TTYPE2 = 'RAWX ' / X-position in Raw Detector Coordinates
TFORM2 = '1B ' / data format of field: BYTE
TUNIT2 = 'pixel ' / physical unit of field
TTYPE3 = 'RAWY ' / Y-position in Raw Detector Coordinates
TFORM3 = '1B ' / data format of field: BYTE
TUNIT3 = 'pixel ' / physical unit of field
TTYPE4 = 'BADFLAG ' / Bad pixel flag
TFORM4 = '16X ' / data format of field: BIT
EXTNAME = 'BADPIX ' / Name of the binary table extension
EXTVER = 4 / There shall be one instance of this extension f
DETNAM = 'DET3 ' / CZT Detector ID (0,1,2 or 3)
TELESCOP= 'NuSTAR ' / Telescope (mission) name
INSTRUME= 'FPMB ' / Instrument name (FPMA or FPMB)
ORIGIN = 'Caltech ' / Source of FITS file
CREATOR = 'FTOOLS 6.10 ' / Creator
VERSION = 1 / Extension version number
FILENAME= 'nuBuserbadpix20100101v001.fits' / File name
CONTENT = 'NuSTAR User Bad Pixel Table' / File content
CCLS0001= 'BCF ' / Daset is a Basic Calibration File
CDTP0001= 'DATA ' / Calibration file contains data
CCNM0001= 'BADPIX ' / Type of calibration data
CVSD0001= '2010-01-01' / Date when this file should first be used
CVST0001= '00:00:00' / Time of day when this file should first be used
CDES0001= 'NuSTAR User Bad Pixel Table' / Description
COMMENT
COMMENT This extension provides, for Detector #3 of the FPMA, the list of
COMMENT user-defined bad pixels.
COMMENT The BADFLAG column contains the bad pixel flags. These are a
COMMENT combination of the following bit settings:
COMMENT
COMMENT b0000000000000001 Bad pixel from on-ground CALDB Bad Pixel File
COMMENT b0000000000000010 Disabled pixel from on-board software
COMMENT b0000000000000100 Bad pixels in the file provided by the user
COMMENT
END

View File

@ -0,0 +1,11 @@
SIMPLE = T / file does conform to FITS standard
BITPIX = 16 / number of bits per data pixel
NAXIS = 0 / number of data axes
EXTEND = T / FITS dataset may contain extensions
COMMENT FITS (Flexible Image Transport System) format is defined in 'Astronomy
COMMENT and Astrophysics', volume 376, page 359; bibcode: 2001A&A...376..359H
TELESCOP= 'NuSTAR ' / Telescope (mission) name
INSTRUME= 'FPMA ' / Instrument name (FPMA or FPMB)
ORIGIN = 'Caltech ' / Source of FITS file
CREATOR = 'FTOOLS 6.10 ' / Creator
END

View File

@ -2,7 +2,7 @@
import numpy as np
import itertools
from pandas import DataFrame, read_csv
from astropy.table import Table, unique
from astropy.table import Table
from astropy.coordinates import SkyCoord
from astropy import units as u
from multiprocessing import get_context, cpu_count
@ -17,7 +17,10 @@ from warnings import filterwarnings
filterwarnings('ignore')
def get_link_list(folder, sort_list=True):
def get_link_list(folder: str, sort_list: bool = True) -> list[str]:
"""
Returns array of paths to all *_cl.evt files in the directory recursively.
"""
links = glob(f'{folder}\\**\\*_cl.evt', recursive=True)
if sort_list:
sorted_list = sorted(links, key=lambda x: stat(x).st_size)
@ -26,7 +29,10 @@ def get_link_list(folder, sort_list=True):
return np.array(links)
def binary_array(num):
def binary_array(num: int) -> list[list[bool]]:
"""
Returns list of all possible combinations of num of bool values.
"""
variants = [[0, 1], ]*num
out = np.zeros((2**num, num), bool)
for idx, level in enumerate(itertools.product(*variants)):
@ -34,8 +40,12 @@ def binary_array(num):
return out
def create_array(filename, mode='Sky'):
temp = fits.getdata(filename, 1)
def create_array(file_path: str, mode: str = 'Sky') -> list[int]:
"""
Returns a 2d array of counts for given observation file.
Modes 'Sky' and 'Det' return arrays in (X,Y) and DET1 respectively.
"""
temp = fits.getdata(file_path, 1)
if mode == 'Sky':
return np.histogram2d(temp['Y'],
temp['X'],
@ -49,6 +59,10 @@ def create_array(filename, mode='Sky'):
def get_wcs(file):
"""
Returns WCS for given observation.
Note that argument here is an opened fits file, not a path.
"""
header = file[1].header
wcs = WCS({
'CTYPE1': header['TCTYP38'], 'CTYPE2': header['TCTYP39'],
@ -61,7 +75,10 @@ def get_wcs(file):
return wcs
def atrous(level=0, max_size=1001):
def atrous(level: int = 0, max_size: int = 1001) -> list[list[float]]:
"""
Returns a trou kernel with the size 2**level and corresponding shape.
"""
base = 1/256*np.array([
[1, 4, 6, 4, 1],
[4, 16, 24, 16, 4],
@ -78,7 +95,19 @@ def atrous(level=0, max_size=1001):
return output
def gauss(level=0, max_size=1000):
def atrous_sig(level: int = 0) -> float:
# sig_values = [0.8908, 0.20066, 0.08551, 0.04122, 0.02042]
sig_values = [0.8725, 0.1893, 0.0946, 0.0473, 0.0237]
if level < 5:
return sig_values[level]
else:
return sig_values[4]/2**(level-4)
def gauss(level: int = 0, max_size: int = 1000) -> list[list[float]]:
"""
Returns gaussian kernel with sigma = 2**level
"""
size = min(5*2**(level+1)+1, max_size)
sigma = 2**(level)
A = 1/(2*np.pi*sigma**2)**0.5
@ -88,6 +117,9 @@ def gauss(level=0, max_size=1000):
def adjecent(array):
"""
Returns two lists of indices of cells adjecent or diagonal to non-zero cells of given array
"""
grid = np.array([
[1, 1, 1],
[1, 0, 1],
@ -99,11 +131,13 @@ def adjecent(array):
np.logical_not(array.mask))
except AttributeError:
output = np.logical_and(output, np.logical_not(array))
output = np.argwhere(output == True)
return output[:, 0], output[:, 1]
return output
def add_borders(array, middle=True):
"""
Returns border mask for an DET1 observation array
"""
mask = np.zeros(array.shape)
datax, datay = np.any(array > 0, 0), np.any(array > 0, 1)
# Add border masks
@ -119,6 +153,9 @@ def add_borders(array, middle=True):
def fill_poisson(array, size_input=32):
"""
Fills all masked elements of an array with poisson signal with local expected value.
"""
if not (isinstance(array, np.ma.MaskedArray)):
print('No mask found')
return array
@ -136,38 +173,69 @@ def fill_poisson(array, size_input=32):
return output
def mirror(array):
size = array.shape[0]
output = np.tile(array, (3, 3))
output[0:size] = np.flipud(output[0:size])
output[2*size:3*size] = np.flipud(output[2*size:3*size])
output[:, 0:size] = np.fliplr(output[:, 0:size])
output[:, 2*size:3*size] = np.fliplr(output[:, 2*size:3*size])
return output
def count_binning(array, count_per_bin: int = 2):
_array = (array[array.mask == False]
if hasattr(array, 'mask') else array)
_array = (_array.flatten()
if array.ndim != 1 else _array)
bin_sum, bin_size = [], []
_sum, _count = 0, 0
for el in _array:
_sum += el
_count += 1
if _sum >= count_per_bin:
bin_sum.append(_sum)
bin_size.append(_count)
_sum, _count = 0, 0
return np.array(bin_sum), np.array(bin_size)
def cstat(expected, data: list, count_per_bin: int = 2) -> float:
_data = data.flatten()
_data = _data[_data.mask == False]
_expected = expected
c_stat = 0
bin_sum_array, bin_count_array = count_binning(_data, count_per_bin)
bin_exp_array = bin_count_array * _expected
c_stat = 2 * (bin_exp_array - bin_sum_array + bin_sum_array * (np.log(bin_sum_array / bin_exp_array))).mean()
return c_stat
class Observation:
def __init__(self, file_name, E_borders=[3,20]):
"""
Main class, contains information about the observation given.
"""
def __init__(self, file_name, E_borders=[3, 20]):
self.filename = file_name
self.name = file_name[file_name.find('nu'):].replace('_cl.evt','')
self.name = file_name[file_name.find('nu'):].replace('_cl.evt', '')
with fits.open(file_name) as file:
self.obs_id = file[0].header['OBS_ID']
self.ra = file[0].header['RA_NOM']
self.dec = file[0].header['DEC_NOM']
self.time_start = file[0].header['TSTART']
self.header = file[0].header
self.exposure = self.header['EXPOSURE']
self.det = self.header['INSTRUME'][-1]
self.wcs = get_wcs(file)
self.data = np.ma.masked_array(*self.get_data(file, E_borders))
self.hard_mask = add_borders(self.data.data, middle=False)
self.exposure = self.header['EXPOSURE']
def get_coeff(self):
"""
Returns normalalizing coefficients for different chips of the observation detector.
Coefficients are obtained from stacked observations in OCC mode.
"""
coeff = np.array([0.977, 0.861, 1.163, 1.05]) if self.det == 'A' else np.array([1.004, 0.997, 1.025, 0.979])
resized_coeff = (coeff).reshape(2, 2).repeat(180, 0).repeat(180, 1)
return resized_coeff
def get_data(self, file, E_borders=[3, 20]):
"""
Returns masked array with DET1 image data for given energy band.
Mask is created from observations badpix tables and to mask the border and gaps.
"""
PI_min, PI_max = (np.array(E_borders)-1.6)/0.04
data = file[1].data.copy()
idx_mask = (data['STATUS'].sum(1) == 0)
@ -181,22 +249,33 @@ class Observation:
mask = np.logical_or(mask, self.get_bad_pix(file))
return output, mask
def get_bad_pix(self, file):
output = np.zeros((360, 360))
kernel = np.ones((5, 5))
for i in range(4):
current_dir = dirname(__file__)
pixpos_file = fits.getdata(f'{current_dir}\\pixpos\\nu{self.det}pixpos20100101v007.fits',i+1)
bad_pix_file = file[3+i].data.copy()
temp = np.zeros(len(pixpos_file), dtype=bool)
for x, y in zip(bad_pix_file['rawx'], bad_pix_file['rawy']):
temp = np.logical_or(temp, np.equal(pixpos_file['rawx'], x)*np.equal(pixpos_file['rawy'], y))
temp = pixpos_file[temp]
output += np.histogram2d(temp['REF_DET1Y'], temp['REF_DET1X'], 360, [[0, 360],[0, 360]])[0]
output = convolve2d(output, kernel, mode='same') > 0
return output
def get_bad_pix(self, file, threshold=0.9):
"""
Creates a mask for observation based on badpix tables.
"""
current_dir = dirname(__file__)
output = np.ones((360, 360))
for det_id in range(4):
badpix = file[3 + det_id].data
badpix_exp = (badpix['TIME_STOP'] - badpix['TIME'])/self.exposure
pixpos = np.load(f'{current_dir}\\pixpos\\ref_pix{self.det}{det_id}.npy', allow_pickle=True).item()
for raw_x, raw_y, exp in zip(badpix['RAWX'], badpix['RAWY'], badpix_exp):
y, x = pixpos[(raw_x, raw_y)]
output[x-3:x+11, y-3:y+11] -= exp
output = np.clip(output, a_min=0, a_max=None)
self.norm_exp_map = output
return output < threshold
def exposure_corr(self, array):
corr = 1 - self.norm_exp_map
corr[corr > 0.1] = 0.
correction_poiss = np.random.poisson(corr*array, corr.shape)
return array + correction_poiss
def wavdecomp(self, mode='gauss', thresh=False, occ_coeff=False):
"""
Performs a wavelet decomposition of image.
"""
# THRESHOLD
if type(thresh) is int:
thresh_max, thresh_add = thresh, thresh/2
@ -205,66 +284,89 @@ class Observation:
# INIT NEEDED VARIABLES
wavelet = globals()[mode]
max_level = 8
conv_out = np.zeros((max_level+1, self.data.shape[0], self.data.shape[1]))
conv_out = np.zeros(
(max_level+1, self.data.shape[0], self.data.shape[1])
)
size = self.data.shape[0]
# PREPARE ORIGINAL DATA FOR ANALYSIS: FILL THE HOLES + MIRROR + DETECTOR CORRECTION
data = self.exposure_corr(self.data)
data = fill_poisson(self.data)
if occ_coeff:
data = data*self.get_coeff()
data = mirror(data)
data_bkg = data.copy()
data = np.pad(data, data.shape[0], mode='reflect')
# ITERATIVELY CONDUCT WAVLET DECOMPOSITION
for i in range(max_level):
conv = fftconvolve(data, wavelet(i), mode='same')
conv[conv < 0] = 0
temp_out = data-conv
# ERRORMAP CALCULATION
if thresh_max != 0:
sig = ((wavelet(i)**2).sum())**0.5
bkg = fftconvolve(data_bkg, wavelet(i), mode='same')
if mode == 'gauss':
sig = ((wavelet(i)**2).sum())**0.5
if mode == 'atrous':
sig = atrous_sig(i)
bkg = fftconvolve(data, wavelet(i), mode='same')
bkg[bkg < 0] = 0
err = (1+np.sqrt(bkg+0.75))*sig
significant = (temp_out > thresh_max*err)[size:2*size, size:2*size]
if thresh_add != 0:
add_significant = (temp_out > thresh_add*err)[size:2*size, size:2*size]
adj = adjecent(significant)
add_condition = np.logical_and(add_significant[adj[0], adj[1]],
np.logical_not(significant[adj[0], adj[1]]))
add_condition = np.logical_and(adj, add_significant)
while (add_condition).any():
to_add = adj[0][add_condition], adj[1][add_condition]
significant[to_add[0], to_add[1]] = True
significant = np.logical_or(significant, add_condition)
adj = adjecent(significant)
add_condition = np.logical_and(add_significant[adj[0], adj[1]],
np.logical_not(significant[adj[0],adj[1]]))
temp_out[size:2*size, size:2*size][np.logical_not(significant)] = 0
add_condition = np.logical_and(adj, add_significant)
significant = np.pad(significant, significant.shape[0], mode='reflect')
temp_out[np.logical_not(significant)] = 0
# WRITING THE WAVELET DECOMP LAYER
conv_out[i] = +temp_out[size:2*size, size:2*size]
# DISCARDING NEGATIVE COMPONENTS OF WAVELETS TO MAKE MASK BY SUMMING WAVELET LAYERS
conv_out[i][conv_out[i] < 0] = 0
data = conv
conv_out[max_level] = conv[size:2*size, size:2*size]
return conv_out
def region_to_raw(self, region):
"""
Returns a hdu_list with positions of masked pixels in RAW coordinates.
"""
x_region, y_region = np.where(region)
tables = []
hdus = []
for i in range(4):
current_dir = dirname(__file__)
pixpos = Table(fits.getdata(f'{current_dir}\\pixpos\\nu{self.det}pixpos20100101v007.fits', i+1))
pixpos = pixpos[pixpos['REF_DET1X'] != -1]
test = np.zeros(len(pixpos['REF_DET1X']), dtype=bool)
for idx, (x, y) in enumerate(zip(pixpos['REF_DET1X'], pixpos['REF_DET1Y'])):
test[idx] = np.logical_and(np.equal(x, x_region), np.equal(y, y_region)).any()
table = Table({'RAWX': pixpos['RAWX'][test], 'RAWY': pixpos['RAWY'][test]})
if not table:
tables.append(table)
else:
tables.append(unique(table))
positions = np.array((pixpos['RAWX'][test], pixpos['RAWY'][test]))
if sum(test) != 0:
positions = np.unique(positions, axis=1)
rawx, rawy = positions[0], positions[1]
time_start = float(78913712)
bad_flag = np.zeros(16, dtype=bool)
bad_flag[13] = 1
columns = []
columns.append(fits.Column('TIME', '1D', 's', array=len(rawx) * [time_start]))
columns.append(fits.Column('RAWX', '1B', 'pixel', array=rawx))
columns.append(fits.Column('RAWY', '1B', 'pixel', array=rawy))
columns.append(fits.Column('BADFLAG', '16X', array=len(rawx) * [bad_flag]))
hdu = fits.BinTableHDU.from_columns(columns)
naxis1, naxis2 = hdu.header['NAXIS1'], hdu.header['NAXIS2']
hdu.header = fits.Header.fromtextfile(f'{current_dir}\\badpix_headers\\nu{self.det}userbadpixDET{i}.txt')
hdu.header['NAXIS1'] = naxis1
hdu.header['NAXIS2'] = naxis2
hdus.append(hdu)
primary_hdu = fits.PrimaryHDU()
primary_hdu.header = fits.Header.fromtextfile(f'{current_dir}\\badpix_headers\\nu{self.det}userbadpix_main.txt')
hdu_list = fits.HDUList([
fits.PrimaryHDU(),
fits.table_to_hdu(tables[0]),
fits.table_to_hdu(tables[1]),
fits.table_to_hdu(tables[2]),
fits.table_to_hdu(tables[3]),
primary_hdu,
*hdus
])
return hdu_list
@ -288,23 +390,32 @@ def process(args):
good_lvl = np.zeros(bin_num, dtype=bool)
good_idx = 0
if obs.exposure > 1000:
wav_obs = obs.wavdecomp('gauss', thresh, occ_coeff=True)
wav_obs = obs.wavdecomp('atrous', thresh, occ_coeff=True)
occ_coeff = obs.get_coeff()
for idx, lvl in enumerate(binary_array(bin_num)):
binary_arr = binary_array(bin_num)
for idx, lvl in enumerate(binary_arr):
try:
region = wav_obs[2:-1][lvl].sum(0) > 0
except ValueError:
region = np.zeros(obs.data.shape, dtype=bool)
masked_obs = np.ma.masked_array(obs.data, mask=region)*occ_coeff
rem_region = np.logical_and(region, np.logical_not(obs.data.mask))
rem_signal[idx] = 1-obs.data[region].sum()/obs.data.sum()
rem_area[idx] = 1 - rem_region.sum()/np.logical_not(obs.data.mask).sum()
poiss_comp[idx] = np.mean((masked_obs-masked_obs.mean())**2/masked_obs.mean())
poiss_comp[idx] = cstat(masked_obs.mean(), masked_obs)
rms[idx] = np.sqrt(((masked_obs-masked_obs.mean())**2).mean())
parameter = lambda idx: ((poiss_comp[idx])**2+((1-rem_area[idx])*0.5)**2)
if (parameter(idx) < parameter(good_idx)):
for idx in range(len(poiss_comp)):
if ((poiss_comp[idx] < poiss_comp[good_idx]) and
(poiss_comp[idx] < poiss_comp[-1] + 0.05) and
(rem_area[idx] > rem_area[-1])):
good_idx = idx
good_lvl = lvl
if good_idx == 0:
good_idx = len(binary_arr) - 1
good_lvl = binary_arr[good_idx]
try:
region = wav_obs[2:-1][good_lvl].sum(0) > 0
if region.sum() > 0:
@ -349,6 +460,12 @@ def process(args):
def process_folder(input_folder=None, start_new_file=None, fits_folder=None, thresh=None):
"""
Generates a fits-table of parameters, folder with mask images in DET1 and BADPIX tables in RAW for all observations in given folder.
Note that observations with exposure < 1000 sec a skipped.
start_new_file can be either 'y' or 'n'.
thresh must be a tuple, e.g. (5,2).
"""
# DIALOGUE
if not (input_folder):
print('Enter path to the input folder')
@ -386,21 +503,44 @@ def process_folder(input_folder=None, start_new_file=None, fits_folder=None, thr
table = {
'obs_id': [], 'detector': [], 'ra': [], 'dec': [],
'lon': [], 'lat': [], 't_start': [], 'exposure': [],
'count_rate': [], 'remaining_area': [], 'poisson_chi2': [],
'poisson_chi2_full': [], 'rms': []
'count_rate': [], 'remaining_area': [], 'poisson_stat': [],
'poisson_stat_full': [], 'rms': []
}
if start_new:
out_table = DataFrame(table)
out_table.to_csv(f'{fits_folder}\\test.csv')
out_table.to_csv(f'{fits_folder}\\test_skipped.csv')
# FILTERING OUT PROCESSED OBSERVATIONS
already_processed_list = read_csv(f'{fits_folder}\\test.csv', index_col=0, dtype={'obs_id':str})
already_skipped_list = read_csv(f'{fits_folder}\\test_skipped.csv', index_col=0, dtype={'obs_id':str})
already_processed = (already_processed_list['obs_id'].astype(str)+already_processed_list['detector']).values
already_skipped = (already_skipped_list['obs_id'].astype(str)+already_skipped_list['detector']).values
obs_list_names = [curr[curr.index('nu')+2:curr.index('_cl.evt')-2] for curr in obs_list]
not_processed = np.array([(curr not in already_processed) for curr in obs_list_names])
not_skipped = np.array([(curr not in already_skipped) for curr in obs_list_names])
already_processed_list = read_csv(
f'{fits_folder}\\test.csv',
index_col=0,
dtype={'obs_id': str}
)
already_skipped_list = read_csv(
f'{fits_folder}\\test_skipped.csv',
index_col=0,
dtype={'obs_id': str}
)
already_processed = (
already_processed_list['obs_id'].astype(str) +
already_processed_list['detector']
).values
already_skipped = (
already_skipped_list['obs_id'].astype(str) +
already_skipped_list['detector']
).values
obs_list_names = [
curr[curr.index('nu')+2:curr.index('_cl.evt')-2]
for curr in obs_list
]
not_processed = np.array([
(curr not in already_processed)
for curr in obs_list_names
])
not_skipped = np.array([
(curr not in already_skipped)
for curr in obs_list_names
])
obs_list = obs_list[np.logical_and(not_processed, not_skipped)]
print(f'Removed already processed observations. {len(obs_list)} observations remain.')
# START PROCESSING
@ -409,8 +549,11 @@ def process_folder(input_folder=None, start_new_file=None, fits_folder=None, thr
for group_idx in range(len(obs_list)//group_size+1):
print(f'Started group {group_idx}')
group_list = obs_list[group_size*group_idx:min(group_size*(group_idx+1), len(obs_list))]
max_size = np.array([stat(file).st_size/2**20 for file in group_list]).max()
process_num = cpu_count() if max_size < 50 else (cpu_count()//2 if max_size < 200 else (cpu_count()//4 if max_size < 1000 else 1))
max_size = np.array([
stat(file).st_size/2**20
for file in group_list
]).max()
process_num = (cpu_count() if max_size < 50 else (cpu_count()//2 if max_size < 200 else (cpu_count()//4 if max_size < 1000 else 1)))
print(f"Max file size in group is {max_size:.2f}Mb, create {process_num} processes")
with get_context('spawn').Pool(processes=process_num) as pool:
packed_args = map(lambda _: (_, thresh), group_list)

Binary file not shown.

Binary file not shown.

Binary file not shown.

Binary file not shown.

Binary file not shown.

Binary file not shown.

Binary file not shown.

Binary file not shown.

View File

@ -5,7 +5,7 @@ with open("README.md", "r") as fh:
setuptools.setup(
name="nuwavsource",
version="0.0.4",
version="0.0.8",
author="Andrey Mukhin",
author_email="amukhin@phystech.edu",
description="A package for source exclusion in NuStar observation data using wavelet decomposition",