Changeset 33

Timestamp:

May 20, 2014 3:23:45 PM (10 years ago)

Author:

JohnLightfoot

Message:

various improvements

File:

• trunk/skygenerator.py (modified) (8 diffs)

trunk/skygenerator.py

@@ -5 +5 @@
 import numpy as np
 import scipy.constants as sc
+
+import common.commonobjects as co
 
 def black_body(temperature, wavenumber):
@@ -38 +40 @@
         spatial_axis = beamsgenerator['spatial axis [arcsec]']
 
-        # oversampling factor
-#        oversample = 10
-        oversample = 1
-
-#        print 'creating sky model with spatial dims [%s,%s] spectral dim %s',\
-#          (oversample * npix, oversample * npix, nspec)
-#        print 'spatial oversampling: %s' % oversample
-        skymodel = np.zeros([npix, npix, nspec], np.float)
+        # skymodel is complex so that its fft can hold truncated version
+        # of infinitesimally sampled map - does that make sense?
+        skymodel = np.zeros([npix, npix, nspec], np.complex)
 
         sky = self.parameters['substages']['Sky']
@@ -83 +80 @@
 
             if type.upper().strip() == 'POINT':
-                print 'creating point source'
                 source_spectrum = self._create_point_source(xpos, ypos,
                   temperature, 
@@ -94 +90 @@
             self.result['sources'][sourcenum]['spectrum'] = source_spectrum
 
-        print 'sources', self.result['sources']
         self.result['sky model'] = skymodel
         self.result['spatial axis'] = spatial_axis
@@ -115 +110 @@
             return
 
-        # calculate fourier phase shift to move point at [0,0] to [rowpos, colpos]
+        # calculate fourier phase shift to move point at [0,0] to 
+        # [rowpos, colpos]
         shiftx = np.zeros([nx], np.complex)
         shiftx[:nx/2] = np.arange(nx/2, dtype=np.complex)
@@ -135 +131 @@
         # sensitivity and source emissivity
         bbmod = np.zeros(np.shape(frequency_axis))
-        print 'frequency axis', frequency_axis
         # ignore floating point errors
         old_settings = np.seterr(all='ignore')
@@ -141 +136 @@
             # simulate real-life 'rounded' cutoffs numerically  
             f1 = 1.0 / (1.0 + pow(cutoffmin/wn, 18) + pow(wn/cutoffmax, 24))
-            print pow(cutoffmin/wn,18), pow(wn/cutoffmax,24)
-            print temperature, wn, black_body(temperature, wn), f1, emissivity
             bbmod[iwn] = black_body(temperature, wn) * f1 * emissivity
         # restore fp behaviour
@@ -163 +156 @@
             skymodel[:,:,iwn] += temp
 
-        return bbmod
+        # return spectrum
+        axis = co.Axis(data=frequency_axis, title='wavenumber', units='cm-1')
+        spectrum = co.Spectrum(data=bbmod, axis=axis, title='Source spectrum',
+          units='W sr-1 m-2 Hz-1')
+
+        return spectrum
 
     def __repr__(self):