Allows frequency-dependent atmosphere

smrt/atmosphere/simple_isotropic_atmosphere.py

@@ -1,5 +1,21 @@
 # coding: utf-8
-"""Implement an isotropic atmosphere with prescribed emission (up and down) and transmittivity
+"""Implement an isotropic atmosphere with prescribed frequency-dependent emission (up and down) and transmittivity.
+
+ TB and transmissivity can be specified as a constant, or a frequency-dependent dictionary
+
+ To make an atmosphere, it is recommended to use the helper function :py:func:`~smrt.atmosphere.simple_isotropic_atmosphere.make_atmosphere`.
+
+
+Examples::
+
+    # the full path import is required
+    from smrt.atmosphere.simple_isotropic_atmosphere import make_atmosphere
+
+    # Constant
+    atmos = make_atmosphere(tbdown=20., tbup=6., trans=1)
+
+    # Frequency-dependent
+    atmos = make_atmosphere(tbdown={10e9: 15.2, 21e9: 23.5})
 
 """
 
@@ -10,16 +26,24 @@
 import numpy as np
 
 # local import
-#from ..core.error import SMRTError
 
 
 #
 # For the developers:
-# A valid Atmosphere is any object that provide a tbdown, tbup and trans functions function of frequency and costheta.
+# A valid atmosphere is any object that provide a tbdown, tbup and trans functions function of frequency and costheta.
 # costheta can be an array or a value, both must be handled.
-# The implementation below is independent of the frequency and the propagation angle
+# The implementation below is independent of the propagation angle
 #
 
+def make_atmosphere(tbdown=0, tbup=0, trans=1):
+
+    """ Construct an atmosphere instance.
+
+    """
+
+    # create the instance
+    return SimpleIsotropicAtmosphere(tbdown=tbdown, tbup=tbup, trans=trans)
+
 
 class SimpleIsotropicAtmosphere(object):
 
@@ -30,11 +54,19 @@ def __init__(self, tbdown=0, tbup=0, trans=1):
         self.constant_trans = trans
 
     def tbdown(self, frequency, costheta, npol):
-        return np.full(len(costheta)*npol, self.constant_tbdown, dtype=np.float32)
+        if isinstance(self.constant_tbdown, dict):
+            return np.full(len(costheta) * npol, self.constant_tbdown[frequency], dtype=np.float32)
+        else:
+            return np.full(len(costheta) * npol, self.constant_tbdown, dtype=np.float32)
 
     def tbup(self, frequency, costheta, npol):
-        return np.full(len(costheta)*npol, self.constant_tbup, dtype=np.float32)
+        if isinstance(self.constant_tbup, dict):
+            return np.full(len(costheta) * npol, self.constant_tbup[frequency], dtype=np.float32)
+        else:
+            return np.full(len(costheta) * npol, self.constant_tbup, dtype=np.float32)
 
     def trans(self, frequency, costheta, npol):
-        return np.full(len(costheta)*npol, self.constant_trans, dtype=np.float32)
-
+        if isinstance(self.constant_trans, dict):
+            return np.full(len(costheta) * npol, self.constant_trans[frequency], dtype=np.float32)
+        else:
+            return np.full(len(costheta) * npol, self.constant_trans, dtype=np.float32)

smrt/atmosphere/test_atmosphere.py

@@ -4,7 +4,7 @@
 
 from smrt import make_snowpack, make_model, sensor, make_soil
 from smrt import sensor_list
-from smrt.atmosphere.simple_isotropic_atmosphere import SimpleIsotropicAtmosphere
+from smrt.atmosphere.simple_isotropic_atmosphere import SimpleIsotropicAtmosphere, make_atmosphere
 
 
 def test_simple_isotropic_atmosphere():
@@ -36,6 +36,12 @@ def test_simple_isotropic_atmosphere():
     assert abs(res1.TbV() - 227.61318467710458) < 1e-2
     assert abs(res2.TbV() - 214.66092232541834) < 1e-2
 
-    #original absorption (Maetzler 1998)
-    #assert abs(res1.TbV() - 223.925496277253) < 1e-2
-    #assert abs(res2.TbV() - 211.7175839349701) < 1e-2
+
+def test_frequency_dependent_atmosphere():
+
+    mu = np.cos(np.arange(0, 90))
+    atmos = make_atmosphere(tbdown={10e9: 15, 21e9: 23}, tbup={10e9:5, 21e9:6}, trans={10e9:1, 21e9:0.95})
+
+    assert np.all(atmos.tbup(frequency=10e9, costheta=mu, npol=2) == 5)
+    assert np.all(atmos.tbdown(frequency=21e9, costheta=mu, npol=2) == 23)
+    assert np.all(atmos.trans(frequency=21e9, costheta=mu, npol=2) == 0.95)