Merge pull request #94 from lopezvoliver/mie

Adding features to re-use MIE subroutine calculations

Py6S/Params/aeroprofile.py

@@ -49,6 +49,24 @@ def PredefinedType(cls, type):
         """
         return "%d" % type
 
+    @classmethod
+    def FromMieFile(cls, mie_file):
+        """Set 6S to use pre-computed outputs of the MIE subroutine written to a text file. 
+           This corresponds to using option iaer=12 in sixs. 
+
+        Arguments:
+
+        * ``mie_file`` -- the name of the text file containing the outputs of the MIE subroutine. 
+
+        Example usage::
+
+          s.aero_profile = AeroProfile.FromMieFile(mie_file = 'FILE.mie')
+
+        """
+        if not mie_file.endswith('.mie'):
+            mie_file = mie_file + '.mie'
+        return "12 from pre-computed mie file\n%s" % mie_file
+
     @classmethod
     def User(cls, **kwargs):
         """Set 6S to use a user-defined aerosol profile based on proportions of standard aerosol components.
@@ -213,7 +231,7 @@ def SunPhotometerDistribution(cls, r, dvdlogr, refr_real, refr_imag):
         comp += " ".join(real) + "\n"
         comp += " ".join(imag) + "\n"
 
-        return header + num + ds + comp + "0 no results saved"
+        return header + num + ds + comp #+ "0 no results saved" 
 
     class AerosolDistribution:
 
@@ -300,7 +318,7 @@ def __str__(self):
                 len(self.values),
             )
             components = "".join(self.values)
-            return result + components + "0 no results saved"
+            return result + components #+ "0 no results saved"
 
     class UserProfile:
 

Py6S/sixs.py

@@ -108,6 +108,10 @@ class SixS(object):
       Note: only one of ``visibility`` or ``aot550`` can be set.
       When setting one, ensure the other is set to ``None``. (By default, ``s.visibility`` is set to None,
       so to set a visibility value you must also do ``s.aot550 = None``)
+
+    * ``mie`` -- Optional name for a text file to save the results of the MIE subroutine. This only applies if 
+    ``s.aero_profile`` is not a ``AeroProfile.User`` or ``AeroProfile.PredefinedType``. (By default, ``s.mie`` 
+      is set to None) 
     """
 
     # Stores the outputs from 6S as an instance of the Outputs class
@@ -155,6 +159,8 @@ def __init__(self, path=None):
 
         self.atmos_corr = AtmosCorr.NoAtmosCorr()
 
+        self.mie = None
+
     def _find_path(self, path=None):
         """Finds the path of the 6S executable.
 
@@ -202,9 +208,33 @@ def _create_atmos_aero_lines(self):
     def _create_aot_vis_lines(self):
         """Create the AOT or Visibility lines for the input file"""
         if not isinstance(self.aero_profile, AeroProfile.UserProfile):
+            # Option to save the output of the MIE subroutine:
+            # only applicable for:
+            # (1) AeroProfile.AerosolDistribution, i.e. it is an instance
+            # of AeroProfile.AerosolDistribution
+            # (2) SunPhotometerDistribution, i.e. it is a string
+            # with more than 2 lines
+            # Not applicable for: 
+            # - PredefinedType (string with only one line)
+            # - User (string with exactly two lines)
+            results_str = "0" 
+            if not isinstance(self.aero_profile, AeroProfile.AerosolDistribution):
+                if isinstance(self.aero_profile, int):
+                    nlines=1
+                else:
+                    nlines = self.aero_profile.count('\n') 
+            else:
+                nlines = 2
+
+            if(nlines>1):
+                if self.mie is not None:
+                    results_str = "1 results saved in %s.mie\n%s\n0" % (self.mie, self.mie)
+                else:
+                    results_str = "0 no results saved\n0"
+
             # We don't need to set AOT or visibility for a UserProfile, but we do for all others
             if self.aot550 is not None:
-                return "0\n%f value\n" % self.aot550
+                return "%s\n%f value\n" % (results_str, self.aot550)
             elif self.visibility is not None:
                 return "%f\n" % self.visibility
             else:

tests/test_profiles.py

@@ -15,6 +15,7 @@
 # You should have received a copy of the GNU Lesser General Public License
 # along with Py6S.  If not, see <http://www.gnu.org/licenses/>.
 
+import os.path
 import unittest
 
 from Py6S import AeroProfile, AtmosProfile, ParameterError, SixS
@@ -74,6 +75,127 @@ def test_aero_profile(self):
                 delta=0.002,
             )
 
+    def test_aero_profile_save_mie(self):
+        s = SixS()
+        s.aero_profile = AeroProfile.SunPhotometerDistribution(
+            [
+                0.050000001,
+                0.065604001,
+                0.086076997,
+                0.112939,
+                0.148184001,
+                0.194428995,
+                0.255104989,
+                0.334715992,
+                0.439173013,
+                0.576227009,
+                0.756052017,
+                0.99199599,
+                1.30157101,
+                1.707757,
+                2.24070191,
+                2.93996596,
+                3.85745192,
+                5.06126022,
+                6.64074516,
+                8.71314526,
+                11.4322901,
+                15,
+            ],
+            [
+                0.001338098,
+                0.007492487,
+                0.026454749,
+                0.058904506,
+                0.082712278,
+                0.073251031,
+                0.040950641,
+                0.014576218,
+                0.003672085,
+                0.001576356,
+                0.002422644,
+                0.004472982,
+                0.007452302,
+                0.011037065,
+                0.014523974,
+                0.016981738,
+                0.017641816,
+                0.016284294,
+                0.01335547,
+                0.009732267,
+                0.006301342,
+                0.003625077,
+            ],
+            [1.47] * 20,
+            [0.0093] * 20,
+        )
+        s.mie='test_mie_file'
+        s.run()
+        self.assertEqual(os.path.exists("test_mie_file.mie"), True)
+
+
+    def test_aero_profile_run_mie(self):
+        s = SixS()
+        s.aero_profile = AeroProfile.SunPhotometerDistribution(
+            [
+                0.050000001,
+                0.065604001,
+                0.086076997,
+                0.112939,
+                0.148184001,
+                0.194428995,
+                0.255104989,
+                0.334715992,
+                0.439173013,
+                0.576227009,
+                0.756052017,
+                0.99199599,
+                1.30157101,
+                1.707757,
+                2.24070191,
+                2.93996596,
+                3.85745192,
+                5.06126022,
+                6.64074516,
+                8.71314526,
+                11.4322901,
+                15,
+            ],
+            [
+                0.001338098,
+                0.007492487,
+                0.026454749,
+                0.058904506,
+                0.082712278,
+                0.073251031,
+                0.040950641,
+                0.014576218,
+                0.003672085,
+                0.001576356,
+                0.002422644,
+                0.004472982,
+                0.007452302,
+                0.011037065,
+                0.014523974,
+                0.016981738,
+                0.017641816,
+                0.016284294,
+                0.01335547,
+                0.009732267,
+                0.006301342,
+                0.003625077,
+            ],
+            [1.47] * 20,
+            [0.0093] * 20,
+        )
+        s.mie='test_mie_file'
+        s.run()
+        app_rad = s.outputs.apparent_radiance
+
+        s.aero_profile = AeroProfile.FromMieFile('test_mie_file')
+        s.run()
+        self.assertAlmostEqual(s.outputs.apparent_radiance, app_rad, delta=0.005)
+
     def test_aero_profile_errors(self):
         with self.assertRaises(ParameterError):
             AeroProfile.User(dust=0.8, oceanic=0.4)