Fixes on reflectance to radiance conversion

README.md

@@ -0,0 +1,15 @@
+sen2cor: Sentinel 2 Level 2A Prototype Processor 
+=================================================================================
+
+sen2cor is a Prototype Processor for processing Sentinel-2 Top of Atmosphere
+reflectance (Level 1C) data into Bottom of Atmoshperic corrected (Level 2A) data.
+It additionally performs a Scene Classification of the corresponding input.
+
+sen2cor works under the Umbrella of Anaconda (Python 2.7)
+on the three platforms Mac OSX, Linux, and Windows 64-bit. This is mandatory! 
+
+This Software contains an IPR of DLR. This part is no Open Source and the
+module for the Atmospheric Correction is attached to this software in form
+of a binary library only.
+
+Please read the docs located at: http://s2tbx.telespazio-vega.de/sen2cor for further info.

sen2cor.egg-info/SOURCES.txt

@@ -4,7 +4,6 @@ README.txt
 SETUP.pdf
 setup.cfg
 setup.py
-sen2cor/L2A_AtmCorr.so
 sen2cor/L2A_Borg.py
 sen2cor/L2A_Config.py
 sen2cor/L2A_Library.py
@@ -146,6 +145,7 @@ sen2cor/cfg/DICO/1.0/SY/representation/representation.xsd
 sen2cor/cfg/DICO/1.0/SY/spatio/spatio.xsd
 sen2cor/cfg/DICO/1.0/SY/tile/tile.xsd
 sen2cor/lib/GlobalSnowMap.tiff
+sen2cor/lib/Readme.txt
 sen2cor/lib/Thumbs.db
 sen2cor/lib/__init__.py
 sen2cor/lib/10/h99000_wv04_rura.atm

sen2cor/L2A_Config.py

@@ -19,8 +19,8 @@ class L2A_Config(Borg):
     _shared = {}
     def __init__(self, sourceDir = False):
         self._processorName = 'Sentinel-2 Level 2A Prototype Processor (Sen2Cor)'
-        self._processorVersion = '2.0.5'
-        self._processorDate = '2015.11.01'
+        self._processorVersion = '2.0.6'
+        self._processorDate = '2015.12.02'
         self._productVersion = '13'
 
         if(sourceDir):
@@ -66,7 +66,6 @@ def __init__(self, sourceDir = False):
             self._adj_km = 1.0
             self._ch940 = array([8,8,9,9,0,0])
             self._cellsize = 0 # pixelsize (m), cellsize (km)
-            self._d2 = None
             self._dem_unit = 0 # [meter] is default DEM heigh unit unit
             self._ibrdf = 0         # brdf correction
             self._thr_g = 0.25      # lower bound for brdf correction
@@ -120,11 +119,11 @@ def __init__(self, sourceDir = False):
             self._wl1900a = array([1.780, 1.970])    # range for interpolation
             self._wv_thr_cirrus = 0.60
             self._sc_lp_blu = 1.0
+            self._d2 = 1.0 # this is a constant, see explanation in line 2061
             self._timestamp = datetime.now()
             self._c0 = None
             self._c1 = None
             self._e0 = None
-            self._d2 = None
             self._wvlsen = None
             self._fwhm = None
             self._acOnly = False
@@ -2059,9 +2058,12 @@ def createL2A_UserProduct(self):
         xp = L2A_XmlParser(self, 'UP1C') 
         pic = xp.getTree('General_Info', 'Product_Image_Characteristics')        
         self._dnScale = float32(pic.QUANTIFICATION_VALUE.text)
-        rc = pic.Reflectance_Conversion
+        # rc = pic.Reflectance_Conversion
         # The earth sun distance correction factor, already squared:
-        self._d2 = 1.0 / float32(rc.U.text)             
+        # Attention! L1C is already corrected by the Earth-Sun distance.
+        # ATCOR expects the uncorrected TOA rad as input and thus applies d2
+        # in order to rescale TOA rad. As this is wrong, d2 is permanently fixed to 1.0 in the constructor
+        # self._d2 = 1.0 / float32(rc.U.text)
 
         if firstInit == True:
             # copy L2A schemes from config_dir into rep_info:    
@@ -2334,18 +2336,6 @@ def timestamp(self, procedure):
         return
 
 
-#     def calcEarthSunDistance2(self, tile):
-#         year =  int(tile[25:29])
-#         month = int(tile[29:31])
-#         day = int(tile[31:33])
-#         doy = date(year,month,day)
-#         doy = int(doy.strftime('%j'))
-#         exc_earth = 0.01673 # earth-sun distance in A.U.
-#         dastr = 1.0 + exc_earth * sin(2 * pi * (doy - 93.5) / 365.)
-#         self._d2 = dastr * dastr
-#         return
-
-
     def parNotFound(self, parameter):
         dummy, basename = os.path.split(self._configFn)
         self.tracer.fatal('Configuration parameter %s not found in %s' % (parameter, basename))

sen2cor/L2A_Tables.py

@@ -1479,12 +1479,12 @@ def TOA_refl2rad(self, index, indataArr):
         '''
         nrows = self.config.nrows
         ncols = self.config.ncols
-
-        # The digital number (DN) as float:         
+
         if(self.config.resolution == 10):
             validBand = [None,0,1,2,None,None,None,3,None,None,None,None,None]
         else:
             validBand = [0,1,2,3,4,5,6,None,7,8,9,10,11]
+
         bandIndex = validBand[index]
         if(bandIndex == None):
             self.config.tracer.debug('Wrong band index %02d for selected resolution %02d', index, self.config.resolution)
@@ -1500,10 +1500,17 @@ def TOA_refl2rad(self, index, indataArr):
         c1  = float32(self.config.c1[bandIndex])
         Es  = float32(self.config.e0[bandIndex])
         sc  = float32(1 / (c1 * self.config.dnScale))
-        # adding again the earth sun correction d2 = 1/u:
+
+        # adding again the earth sun correction d2 = 1.0/U:
         pi32_d2 = float32(pi) * self.config.d2
-        rad     = zeros_like(rho)
+
+        # The final formaula is:
+        # rad = rho * cos(radians(sza)) * Es * sc / (pi * d2)
+        # where: d2 = 1.0 / U
+        # scale: 1 / (0.001 * 1000) = 1 (default)
+
         # To reduce the memory consumption for 10m images:
+        rad = zeros_like(rho)
         for i in range(nrows):
             rad_sza = float32(radians(sza[i,:]))
             cos_sza = float32(cos(rad_sza))

setup.py

@@ -22,8 +22,8 @@
 #
 # This needs to be changed with each new version:
 #------------------------------------------------
-version = '2.0.5'
-longVersion = '02.00.05'
+version = '2.0.6'
+longVersion = '02.00.06'
 #
 # Do not change anything below:
 #----------------------------------------------------------------------------