ComboCode.cc.data.instruments.Spire

23 24 """ 25 Environment with several tools for Spire molecular spectroscopy managing. 26 27 """ 28

29 - def __init__(self,star_name,resolution,oversampling,\ 30 path=None,intrinsic=1,path_linefit=''):

31 32 ''' 33 Initializing an instance of Spire(). 34 35 @param star_name: Name of the star from Star.dat 36 @type star_name: string 37 @param resolution: The spectral resolution of the SPIRE apodized 38 spectra 39 @type resolution: float 40 @param oversampling: The SPIRE instrumental oversampling, for correct 41 convolution of the Sphinx output. This is the 42 oversampling factor for the apodized spectra! 43 @type oversampling: int 44 45 @keyword path: Output folder in the code's home folder 46 47 (default: 'codeSep2010') 48 @type path: string 49 @keyword intrinsic: Use the intrinsic Sphinx line profiles for 50 convolving with the spectral resolution? Otherwise 51 the beam convolved line profiles are used. 52 53 (default: 1) 54 @type intrinsic: bool 55 @keyword path_linefit: The folder name for linefit results from Hipe 56 (created by Maries script, assuming her syntax). 57 The folder is located in 58 $path_pacs$/$star_name$/. If no folder is given, 59 no linefits are processed. 60 61 (default: '') 62 @type path_linefit: string 63 64 ''' 65 66 super(Spire,self).__init__(star_name=star_name,code='GASTRoNOoM',\ 67 path=path,oversampling=oversampling,\ 68 path_linefit=path_linefit,blend_factor=1.3,\ 69 instrument_name='SPIRE',intrinsic=intrinsic) 70 #- resolution is given in cm^-1 71 self.resolution = float(resolution) 72 self.sigma = self.resolution/(2.*sqrt(2.*log(2.))) 73 self.sphinx_convolution = dict() 74 if not self.resolution: 75 print 'WARNING! SPIRE resolution is undefined!' 76 self.readLineFit() 77 78 #-- Convenience path 79 cc.path.gout = os.path.join(cc.path.gastronoom,self.path)

80 81

82 - def prepareSphinx(self,star_grid,redo_sphinx_prep=0):

83 84 ''' 85 86 Prepare Sphinx SPIRE models by checking if the convolved spectrum is 87 already present in the Spire object and if not calling convolveSphinx. 88 89 @param star_grid: list of Star() instances 90 @type star_grid: list[Star()] 91 @keyword redo_sphinx_prep: redo the sphinx prep, regardless of whether 92 this Spire instance already did it once, 93 for instance in case the star_grid changed 94 95 (default: 0) 96 @type redo_sphinx_prep: bool 97 98 ''' 99 100 if not self.sphinx_convolution or redo_sphinx_prep: 101 print 'Convolving with Gaussian and rebinning to data ' + \ 102 'wavelength grid.' 103 for i,star in enumerate(star_grid): 104 print '* Sphinx model %i out of %i.' %(i+1,len(star_grid)) 105 if not star['LAST_GASTRONOOM_MODEL']: 106 print '* No cooling model found.' 107 else: 108 self.sphinx_convolution[i] = dict() 109 star['LAST_SPIRE_MODEL'] = i 110 self.__convolveSphinx(star=star) 111 if self.sphinx_convolution[i]: 112 print '* Model %i with cooling id %s is done!'\ 113 %(i,star['LAST_GASTRONOOM_MODEL']) 114 else: 115 star['LAST_SPIRE_MODEL'] = None

116 117 118

119 - def getSphinxConvolution(self,star,fn):

120 121 ''' 122 Return the sphinx convolution and return if it has already been done. 123 124 Returns None if the convolution is not available. 125 126 @param star: The Star() object 127 @type star: Star() 128 @param fn: The filename of the dataset (band) for which the convolution 129 is to be returned. 130 @type fn: str 131 132 @return: The sphinx convolution result. (wavelength,flux) 133 @rtype: array 134 135 ''' 136 137 if star['LAST_SPIRE_MODEL'] is None: 138 return ([],[]) 139 ifn = self.data_filenames.index(fn) 140 return (self.data_wave_list[ifn],\ 141 self.sphinx_convolution[star['LAST_SPIRE_MODEL']][fn])

142 143 144 145

146 - def __convolveSphinx(self,star):

147 148 ''' 149 Convolve the Sphinx output with the SPIRE resolution. The convolution 150 is done in wave number (cm^-1). 151 152 @param star: The Star() object for which Sphinx profiles are loaded 153 @type star: Star() 154 155 ''' 156 157 #- Get sphinx model output and merge, for all star models in star_grid 158 if not self.resolution: 159 print '* Resolution is undefined. Cannot convolve Sphinx.' 160 return 161 print '* Reading Sphinx model and merging.' 162 sphinx_wav,sphinx_flux = star['LAST_GASTRONOOM_MODEL'] \ 163 and self.mergeSphinx(star) \ 164 or [[],[]] 165 if not sphinx_wav: 166 print '* No Sphinx data found.' 167 return 168 sphinx_wav = 1./array(sphinx_wav)*10**(4) 169 sphinx_flux = array(sphinx_flux) 170 sphinx_wav = sphinx_wav[::-1] 171 sphinx_flux = sphinx_flux[::-1] 172 173 #-- eliminate some of the zeroes in the grid to reduce calculation time 174 # (can reduce the array by a factor up to 100!!) 175 s = self.sigma 176 lcs = array(sorted([1./line.wavelength 177 for line in star['GAS_LINES']])) 178 new_wav, new_flux = [sphinx_wav[0]],[sphinx_flux[0]] 179 for w,f in zip(sphinx_wav[1:],sphinx_flux[1:]): 180 if f != 0 or (w < 5*s+lcs[argmin(abs(lcs-w))] \ 181 and w > lcs[argmin(abs(lcs-w))]-5*s): 182 new_wav.append(w) 183 new_flux.append(f) 184 new_wav, new_flux = array(new_wav), array(new_flux) 185 186 #-- convolve the model fluxes with a gaussian and constant sigma(spire) 187 print '* Convolving Sphinx model for SPIRE.' 188 convolution = Data.convolveArray(new_wav,new_flux,s) 189 190 for data_wav,fn in zip(self.data_wave_list,self.data_filenames): 191 rebinned = [] 192 #-- Convert wavelengths to wave number for integration, and reverse 193 data_cm = data_wav[::-1] 194 data_cm = 1./data_cm*10**4 195 rebinned = [trapz(y=convolution[abs(new_wav-wavi)<=self.resolution/self.oversampling],\ 196 x=new_wav[abs(new_wav-wavi)<=self.resolution/self.oversampling])\ 197 /(self.resolution/self.oversampling) 198 for wavi in data_cm] 199 #-- Reverse the rebinned fluxes so they match up with the 200 # wavelength grid. 201 rebinned = array(rebinned)[::-1] 202 self.sphinx_convolution[star['LAST_SPIRE_MODEL']][fn] = rebinned

203 204 205

206 - def readLineFit(self):

207 208 ''' 209 Read the data from the line fit procedure done with Maries Hipe 210 script. 211 212 Assumes structure and syntax as given in the example file 213 /home/robinl/Data/SPIRE/rdor/lineFit/lineFitResults 214 215 The line fit results are saved in self.linefit as a np.recarray. 216 217 The relative FWHM is given with respect to the intrinsic SPIRE 218 resolution for unapodized spectra, i.e. 0.04 cm^-1. 219 220 The columns include (with unit indicated): 221 band 222 freq_in (GHz), 223 freq_fit (GHz), 224 wave_fit (micron), 225 line_flux (W/m2), 226 line_flux_err (W/m2), 227 line_flux_rel (ratio), 228 line_peak (W/m2 m), 229 fwhm_fit_freq (GHz), 230 fwhm_fit (micron), 231 fwhm_rel 232 233 ''' 234 235 dd = super(Spire,self).readLineFit(**dict([('start_row',2)])) 236 if dd is None: return 237 238 #-- Remove +/- 239 del dd[5] 240 names = ('band','freq_in','freq_fit','wave_fit','line_flux',\ 241 'line_flux_err','line_flux_rel','line_peak','fwhm_fit_ghz',\ 242 'fwhm_fit', 'fwhm_rel') 243 self.linefit = np.recarray(shape=(len(dd[-1]),),\ 244 dtype=zip(names,['|S3']+[float]*10)) 245 for n,d in zip(names,dd): 246 self.linefit[n] = d

Source Code for Module ComboCode.cc.data.instruments.Spire