Hide keyboard shortcuts

Hot-keys on this page

r m x p   toggle line displays

j k   next/prev highlighted chunk

0   (zero) top of page

1   (one) first highlighted chunk

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

95

96

97

98

99

100

101

102

103

104

105

106

107

108

109

110

111

112

113

114

115

116

117

118

119

120

121

122

123

124

125

126

127

128

129

130

131

132

133

134

135

136

137

138

139

140

141

142

143

144

145

146

147

148

149

150

151

152

153

154

155

156

157

158

159

160

161

162

163

164

165

166

167

168

169

170

171

172

173

174

175

176

177

178

179

180

181

182

183

184

185

186

187

188

189

190

191

192

193

194

195

196

197

198

199

200

201

202

203

204

205

206

207

208

209

import os 

 

from .masks import MaskHelper 

from .target import Target 

from .utils import astropyHeaderFromDict, inheritDocstrings 

from .wavelengthArray import WavelengthArray 

 

__all__ = ["PfsSimpleSpectrum"] 

 

 

@inheritDocstrings 

class PfsSimpleSpectrum: 

"""Spectrum for a single object 

 

This base class is suitable for model spectra which have not been extracted 

from observations. 

 

Parameters 

---------- 

target : `pfs.datamodel.Target` 

Target information. 

wavelength : `numpy.ndarray` of `float` 

Array of wavelengths. 

flux : `numpy.ndarray` of `float` 

Array of fluxes. 

mask : `numpy.ndarray` of `int` 

Array of mask pixels. 

flags : `pfs.datamodel.MaskHelper` 

Helper for dealing with symbolic names for mask values. 

metadata : `dict` (`str`: POD), optional 

Keyword-value pairs for the header. 

""" 

filenameFormat = None # Subclasses should override 

 

def __init__(self, target, wavelength, flux, mask, flags, metadata=None): 

self.target = target 

self.wavelength = wavelength 

self.flux = flux 

self.mask = mask 

self.flags = flags 

self.metadata = metadata if metadata is not None else {} 

 

self.length = len(wavelength) 

self.validate() 

 

def validate(self): 

"""Validate that all the arrays are of the expected shape""" 

assert self.wavelength.shape == (self.length,) 

assert self.flux.shape == (self.length,) 

assert self.mask.shape == (self.length,) 

 

def __len__(self): 

"""Return the length of the arrays""" 

return self.length 

 

def getIdentity(self): 

"""Return the identity of the spectrum 

 

Returns 

------- 

identity : `dict` 

Key-value pairs that identify this spectrum. 

""" 

return self.target.identity 

 

@classmethod 

def _readImpl(cls, fits): 

"""Implementation for reading from FITS file 

 

Parameters 

---------- 

fits : `astropy.io.fits.HDUList` 

Opened FITS file. 

 

Returns 

------- 

kwargs : ``dict`` 

Keyword arguments for constructing spectrum. 

""" 

data = {} 

data["flux"] = fits["FLUX"].data 

data["mask"] = fits["MASK"].data 

 

# Wavelength can be specified in an explicit extension, or as a WCS in the header 

if "WAVELENGTH" in fits: 

wavelength = fits["WAVELENGTH"].data 

else: 

wavelength = WavelengthArray.fromFitsHeader(fits["FLUX"].header, len(fits["FLUX"].data)) 

data["wavelength"] = wavelength 

 

data["flags"] = MaskHelper.fromFitsHeader(fits["MASK"].header) 

data["target"] = Target.fromFits(fits) 

return data 

 

@classmethod 

def readFits(cls, filename): 

"""Read from FITS file 

 

This API is intended for use by the LSST data butler, which handles 

translating the desired identity into a filename. 

 

Parameters 

---------- 

filename : `str` 

Filename of FITS file. 

 

Returns 

------- 

self : ``cls`` 

Constructed instance, from FITS file. 

""" 

import astropy.io.fits 

with astropy.io.fits.open(filename) as fd: 

data = cls._readImpl(fd) 

return cls(**data) 

 

@classmethod 

def read(cls, identity, dirName="."): 

"""Read from file 

 

This API is intended for use by science users, as it allows selection 

of the correct file from parameters that make sense, such as which 

catId, objId, etc. 

 

Parameters 

---------- 

identity : `dict` 

Keyword-value pairs identifying the data of interest. Common keywords 

include ``catId``, ``tract``, ``patch``, ``objId``. 

dirName : `str`, optional 

Directory from which to read. 

 

Returns 

------- 

self : ``cls`` 

Spectrum read from file. 

""" 

filename = os.path.join(dirName, cls.filenameFormat % identity) 

return cls.readFits(filename) 

 

def _writeImpl(self, fits): 

"""Implementation for writing to FITS file 

 

We attempt to write the wavelength to the header (as a WCS; this results 

in a modest size savings), which works if the wavelength is a specified 

as a `WavelengthArray`; otherwise we write it as an explicit extension. 

 

Parameters 

---------- 

fits : `astropy.io.fits.HDUList` 

List of FITS HDUs. This has a Primary HDU already, the header of 

which may be supplemented with additional keywords. 

 

Returns 

------- 

header : `astropy.io.fits.Header` 

FITS headers which may contain the wavelength WCS. 

""" 

from astropy.io.fits import ImageHDU, Header 

haveWavelengthHeader = False 

try: 

header = self.wavelength.toFitsHeader() # For WavelengthArray 

haveWavelengthHeader = True 

except AttributeError: 

header = Header() 

if self.metadata: 

header.extend(astropyHeaderFromDict(self.metadata)) 

fits.append(ImageHDU(self.flux, header=header, name="FLUX")) 

maskHeader = astropyHeaderFromDict(self.flags.toFitsHeader()) 

maskHeader.extend(header) 

fits.append(ImageHDU(self.mask, header=maskHeader, name="MASK")) 

if not haveWavelengthHeader: 

fits.append(ImageHDU(self.wavelength, header=header, name="WAVELENGTH")) 

self.target.toFits(fits) 

return header 

 

def writeFits(self, filename): 

"""Write to FITS file 

 

This API is intended for use by the LSST data butler, which handles 

translating the desired identity into a filename. 

 

Parameters 

---------- 

filename : `str` 

Filename of FITS file. 

""" 

from astropy.io.fits import HDUList, PrimaryHDU 

fits = HDUList() 

fits.append(PrimaryHDU()) 

self._writeImpl(fits) 

with open(filename, "wb") as fd: 

fits.writeto(fd) 

 

def write(self, dirName="."): 

"""Write to file 

 

This API is intended for use by science users, as it allows setting the 

correct filename from parameters that make sense, such as which 

catId, objId, etc. 

 

Parameters 

---------- 

dirName : `str`, optional 

Directory to which to write. 

""" 

identity = self.getIdentity() 

filename = os.path.join(dirName, self.filenameFormat % identity) 

return self.writeFits(filename)