Merge branch 'dev'

Author: oceancolorcoder

Main.py

@@ -16,10 +16,6 @@
 import sys
 import time
 from PyQt5 import QtCore, QtGui, QtWidgets
-from PIL import Image, ImageTk
-
-# import requests
-# from tqdm import tqdm
 
 from Source import PACKAGE_DIR as CODE_HOME
 from Source.MainConfig import MainConfig
@@ -31,7 +27,7 @@
 from Source.SeaBASSHeaderWindow import SeaBASSHeaderWindow
 from Source.Utilities import Utilities
 
-VERSION = "1.2.9"
+VERSION = "1.2.10"
 
 
 class Window(QtWidgets.QWidget):

README.md

@@ -16,7 +16,7 @@ retrieval. Data outputs are formatted to text files for submission to the [SeaBA
 Currently, HyperCP supports <a href='https://www.seabird.com/'>Sea-Bird Scientific</a> HyperSAS packages with and
 without SolarTracker or pySAS robotic platforms as well as [TriOS](https://www.trios.de/en/radiometers.html) used in manual configuration. If you are interested in integrating support for your platform, contact us at the email addresses below or in the Discussions tab of the GitHub repository.
 
-## Version 1.2.9
+## Version 1.2.10
 
 ```
 The MIT license

Source/ConfigWindow.py

@@ -220,7 +220,7 @@ def initUI(self):
         l1bSublabel6 = QtWidgets.QLabel("    Fallback values when no model available:", self)
         # l1bSublabel5.setOpenExternalLinks(True)
         self.l1bGetAncCheckBox1 = QtWidgets.QCheckBox("GMAO MERRA2", self)
-        self.l1bGetAncCheckBox2 = QtWidgets.QCheckBox("ECMWF EAC4", self)
+        self.l1bGetAncCheckBox2 = QtWidgets.QCheckBox("ECMWF CAMS", self)
 
         # If clicked trigger l1bGetAncCheckBoxUpdate
         self.l1bGetAncCheckBox1.clicked.connect(lambda: self.l1bGetAncCheckBoxUpdate('NASA_Earth_Data'))
@@ -1730,7 +1730,7 @@ def l1bGetAncCheckBoxUpdate(self,ancillarySource):
             credentials.credentialsWindow('NASA_Earth_Data')
             self.l1bGetAncUntickIfNoCredentials('NASA_Earth_Data')
         elif self.l1bGetAncCheckBox2.isChecked():
-            print("ConfigWindow - l1bGetAncCheckBoxUpdate ECMWF EAC4")
+            print("ConfigWindow - l1bGetAncCheckBoxUpdate ECMWF CAMS")
             ConfigFile.settings["bL1bGetAnc"] = 2
             credentials.credentialsWindow('ECMWF_ADS')
             self.l1bGetAncUntickIfNoCredentials('ECMWF_ADS')

Source/GetAnc_credentials.py

@@ -39,7 +39,7 @@ def credentialsSpec(credentialsID):
         out['URL_string'] = 'https://urs.earthdata.nasa.gov'
     elif credentialsID == 'ECMWF_ADS':
         out['credentials_filename'] = '.ecmwf_ads_credentials.json'
-        out['title'] = 'Login ECMWF ADS credentials - EAC-4 Ancillary'
+        out['title'] = 'Login ECMWF ADS credentials - CAMS Ancillary'
         out['key_string'] = 'url'
         out['secret_string'] = 'key'
         out['URL_string'] = 'https://ads.atmosphere.copernicus.eu/how-to-api'

Source/GetAnc_ecmwf.py

@@ -19,34 +19,18 @@
 
 class GetAnc_ecmwf:
 
-    def timeStamp2yrMnthDayHrMinSec(timestamp):
-        '''
-        Definition for timestamp managing
-        :param timeStamp: a string, the time in UTC with format yyyy-mm-ddThh:MM:ss
-        '''
-        date = timestamp.split('T')[0]
-        time = timestamp.split('T')[1]
-
-        year = date.split('-')[0]
-        month = date.split('-')[1]
-        day = date.split('-')[2]
-
-        hour = time.split(':')[0]
-        minute = time.split(':')[1]
-        second = time.split(':')[2]
-
-        return year,month,day,hour,minute,second
-
     def ECMWF_latLonTimeTags(lat, lon, timeStamp, latRes, lonRes, timeResHours):
         '''
         :param timeStamp: a string, the time in UTC with format yyyy-mm-ddThh:MM:ss
         :param lat: a float, the query latitude in degrees North
         :param lon: a float, the query longitude in degrees East
+        :param timeResHours: an integer, time resolution of the queried ECMWF dataset
         :return:
         latEff: a float, the effective latitude - i.e. to the closest resolution degree
         lonEff: a float, the effective latitude - i.e. to the closest resolution degree
         latLonTag: a string, a tag to indicate the effective lat/lon.
-        dateTag: a string, a tag to indicate time stamp, to the hour
+        dateTagEff: a string, a tag to indicate effective date (i.e. date in UTC rounded to closest hour)
+        timeStampEff: a string, a tag to indicate effective time (i.e. time in UTC rounded to closest hour)
         '''
 
         # Lat-Lon
@@ -63,31 +47,19 @@ def ECMWF_latLonTimeTags(lat, lon, timeStamp, latRes, lonRes, timeResHours):
                                     str(int(np.sign(latEff))).replace('-1', 'S').replace('1', 'N'),
                                     np.abs(latEff * (10 ** latSigFigures)))
 
-        # Time; choose 00:00 or 12:00 model
-        year, month, day, hour, minute, second = GetAnc_ecmwf.timeStamp2yrMnthDayHrMinSec(timeStamp)
-        dTtimeStamp = datetime.datetime(int(year), int(month), int(day), int(hour), int(minute), int(second),tzinfo=datetime.timezone.utc)
-        dTtimeVec = [datetime.datetime(int(year), int(month), int(day), int(0), int(0), int(0),tzinfo=datetime.timezone.utc),\
-                     datetime.datetime(int(year), int(month), int(day), int(12), int(0), int(0),tzinfo=datetime.timezone.utc)]
-        timeEff = min(dTtimeVec, key=lambda x: abs(x - dTtimeStamp))
-
-        # timeSec = datetime.datetime(int(year), int(month), int(day), int(hour), int(minute), int(second),tzinfo=datetime.timezone.utc).timestamp()
-        # timeEffSec = np.round(timeSec / (3600 * timeResHours)) * (3600 * timeResHours) - 7200
-        # timeStampEff = str(pd.Timestamp(datetime.datetime.fromtimestamp(timeEffSec))).replace(' ','T')
-        timeStampEff = str(pd.Timestamp(timeEff)).replace(' ','T')
-
-        dateTag = timeStampEff.replace(':','').replace('-','').replace('T','')
-        # print(timeStamp)
-        # print(year, month, day, hour, minute, second)
-        # print(timeSec)
-        # print(timeEffSec)
-        # print(timeStampEff)
+        # Convert to a datetime object
+        epoch_time = datetime.datetime.strptime(':'.join(timeStamp.split(':')[:-2]), '%Y-%m-%dT%H:%M:%S').timestamp()
+        timeResHoursSecs = 3600 * timeResHours
+        rounded_epoch_time = round(epoch_time / timeResHoursSecs) * timeResHoursSecs
+        rounded_timestamp = datetime.datetime.fromtimestamp(rounded_epoch_time).strftime('%Y-%m-%dT%H:%M:%S')
+        dateTagEff, timeStampEff = rounded_timestamp.split('T')
 
-        return latEff,lonEff,timeStampEff,latLonTag,dateTag
+        return latEff,lonEff,latLonTag,dateTagEff,timeStampEff
 
 
-    def EAC4_download_ensembles(lat, lon, timeStamp, EAC4_variables, pathOut):
+    def CAMS_download_ensembles(lat, lon, dateTag, timeTag, CAMS_variables, pathOut):
         '''
-        Performs CDSAPI command to download the required data from EAC4 (dataset "cams-global-atmospheric-composition-forecasts") in netCDF
+        Performs CDSAPI command to download the required data from CAMS (dataset "cams-global-atmospheric-composition-forecasts") in netCDF
         format. It will retrieve the variables in a single space-time point corresponding to
 
         For more information, please check: https://ads.atmosphere.copernicus.eu/cdsapp#!/dataset/cams-global-atmospheric-composition-forecasts?tab=overview
@@ -96,38 +68,42 @@ def EAC4_download_ensembles(lat, lon, timeStamp, EAC4_variables, pathOut):
         :param lat: a float, the query latitude in degrees North
         :param lon: a float, the query longitude in degrees East
 
-        :param EAC4_variables: a list, with the variables of interest
+        :param CAMS_variables: a list, with the variables of interest
         :param pathOut: full path to output, a netCDF file with the requested variables.
         :return:
         '''
 
         if os.path.exists(pathOut):
             pass
         else:
-            print(f'Nearest model found at {timeStamp}')
             url,key = read_user_credentials('ECMWF_ADS')
 
-            year, month, day, hour, _, _ = GetAnc_ecmwf.timeStamp2yrMnthDayHrMinSec(timeStamp)
+            year = dateTag.split('-')[0]
+            hour = timeTag.split(':')[0]
 
-            if int(year) < 2003:
-                print('EAC4 dataset not available before 2003, skipping')
+            hourForecast = '%02d' % (int(hour) // 12)
+            leadtime     = str(int(hour) % 12)
+
+            if int(year) < 2015:
+                print('CAMS dataset not available before 2015, skipping')
             else:
                 try:
                     c = cdsapi.Client(timeout=5, url=url, key=key)
                     c.retrieve(
                         'cams-global-atmospheric-composition-forecasts',
                         {
-                            'format': 'netcdf',
                             'type' : 'forecast',
-                            'variable': list(EAC4_variables.keys()),
-                            'date': '%s-%s-%s/%s-%s-%s' % (year, month, day, year, month, day),
-                            'time': '%s:00' % (hour,),
+                            'variable': list(CAMS_variables.keys()),
+                            'date': '%s/%s' % (dateTag, dateTag),
                             'area': [lat, lon, lat, lon],
-                            'leadtime_hour': '0',
+                            'time': '%s:00' % hourForecast,
+                            'leadtime_hour': leadtime,
+                            'format': 'netcdf',
+                            'download_format': 'unarchived',
                         },
                         pathOut)
                 except:
-                    print('EAC4 atmospheric data could not be retrieved. Check inputs.')
+                    print('CAMS atmospheric data could not be retrieved. Check inputs.')
                     exit()
 
 
@@ -152,49 +128,49 @@ def get_ancillary_main(lat, lon, timeStamp, pathAncillary):
 
         ancillary = {}
 
-        #################### EAC4 ####################
-        pathEAC4 = pathAncillary
-        if not os.path.exists(pathEAC4):
-            os.mkdir(pathEAC4)
+        #################### CAMS ####################
+        pathCAMS = pathAncillary
+        if not os.path.exists(pathCAMS):
+            os.mkdir(pathCAMS)
 
-        EAC4_variables = {
+        CAMS_variables = {
         '10m_u_component_of_wind' :'u10',
         '10m_v_component_of_wind' :'v10',
         'total_aerosol_optical_depth_550nm' :'aod550'
         }
 
-        EAC4nc = {}
+        CAMSnc = {}
 
 
         # Check https://ads.atmosphere.copernicus.eu/cdsapp#!/dataset/cams-global-atmospheric-composition-forecasts?tab=overview
         latRes = 0.4
         lonRes = 0.4
-        timeResHours = 12
+        timeResHours = 1
 
-        latEff, lonEff, timeStampEff, latLonTag, dateTag = GetAnc_ecmwf.ECMWF_latLonTimeTags(lat, lon, timeStamp, latRes, lonRes, timeResHours)
+        latEff, lonEff, latLonTag, dateTagEff, timeStampEff = GetAnc_ecmwf.ECMWF_latLonTimeTags(lat, lon, timeStamp, latRes, lonRes, timeResHours)
 
-        pathOut = os.path.join(pathEAC4, 'EAC4_%s_%s.nc' % (latLonTag, dateTag))
+        pathOut = os.path.join(pathCAMS, 'CAMS_%s_%s_%s.nc' % (latLonTag, dateTagEff.replace('-',''), timeStampEff.replace(':','')))
 
-        GetAnc_ecmwf.EAC4_download_ensembles(latEff, lonEff, timeStampEff, EAC4_variables, pathOut)
+        GetAnc_ecmwf.CAMS_download_ensembles(latEff, lonEff, dateTagEff, timeStampEff, CAMS_variables, pathOut)
 
         try:
-            EAC4nc['reanalysis'] = xr.open_dataset(pathOut,engine='netcdf4')
-            EAC4_flag = True
+            CAMSnc['reanalysis'] = xr.open_dataset(pathOut,engine='netcdf4')
+            CAMS_flag = True
         except:
-            EAC4_flag = False
-            print('EAC4 data missing. Skipping...')
+            CAMS_flag = False
+            print('CAMS data missing. Skipping...')
 
-        if EAC4_flag:
+        if CAMS_flag:
             try:
-                for EAC4_variable, shortName in EAC4_variables.items():
-                    var = EAC4nc['reanalysis'][shortName]
-                    ancillary[EAC4_variable] = {}
-                    ancillary[EAC4_variable]['value']      = var.values[0][0][0]
-                    ancillary[EAC4_variable]['units']      = var.units
-                    ancillary[EAC4_variable]['long_name']  = var.long_name
-                    ancillary[EAC4_variable]['source']     = 'EAC4 (ECMWF). https://ads.atmosphere.copernicus.eu/cdsapp#!/dataset/cams-global-atmospheric-composition-forecasts?tab=overview'
+                for CAMS_variable, shortName in CAMS_variables.items():
+                    var = CAMSnc['reanalysis'][shortName]
+                    ancillary[CAMS_variable] = {}
+                    ancillary[CAMS_variable]['value']      = var.values[0][0][0][0]
+                    ancillary[CAMS_variable]['units']      = var.units
+                    ancillary[CAMS_variable]['long_name']  = var.long_name
+                    ancillary[CAMS_variable]['source']     = 'CAMS (ECMWF). https://ads.atmosphere.copernicus.eu/cdsapp#!/dataset/cams-global-atmospheric-composition-forecasts?tab=overview'
             except:
-                print('Problem processing EAC4 data. Skipping...')
+                print('Problem processing CAMS data. Skipping...')
 
         return ancillary
 
@@ -223,11 +199,11 @@ def getAnc_ecmwf(inputGroup):
             ancillary = GetAnc_ecmwf.get_ancillary_main(lat[index], lon[index], lat_timeStamp, ancPath)
 
             # position retrieval index has been confirmed manually in SeaDAS
-            uWind = ancillary['10m_u_component_of_wind']['value'][0]
-            vWind = ancillary['10m_v_component_of_wind']['value'][0]
+            uWind = ancillary['10m_u_component_of_wind']['value']
+            vWind = ancillary['10m_v_component_of_wind']['value']
             modWind.append(np.sqrt(uWind*uWind + vWind*vWind)) # direction not needed
             #ancAOD = aerGroup.getDataset("TOTEXTTAU")
-            modAOD.append(ancillary['total_aerosol_optical_depth_550nm']['value'][0])
+            modAOD.append(ancillary['total_aerosol_optical_depth_550nm']['value'])
 
         modData = HDFRoot()
         modGroup = modData.addGroup('ECMWF')

Source/ProcessL1aTriOS.py

@@ -360,6 +360,7 @@ def formatting_instrument(name, cal_path, input_file, root, configPath):
         c3 = float(gp.attributes['c3s'])
         wl = []
         for i in range(1,256):
+        # for i in range(1,data.shape[1]+1):
             wl.append(str(round((c0 + c1*(i+1) + c2*(i+1)**2 + c3*(i+1)**3), 2)))
 
         #Create Data (LI,LT,ES) dataset

Source/ProcessL1b.py

@@ -162,7 +162,10 @@ def read_unc_coefficient_frm(root, inpath, classbased_dir):
         Utilities.RenameUncertainties_FullChar(root)
 
         # interpolate LAMP and PANEL to full wavelength range
-        Utilities.interpUncertainties_FullChar(root)
+        success = Utilities.interpUncertainties_FullChar(root)
+        if not success:
+            print('interpUncertainties_FullChar failed.')
+            return None
 
         # generate temperature coefficient
         Utilities.UncTempCorrection(root)

Source/ProcessL1bqc.py

@@ -662,7 +662,7 @@ def QC(node):
                     print(msg)
                     Utilities.writeLogFile(msg)
                     return False
-                
+
                 if py6sGroup is not None:
                     Utilities.filterData(py6sGroup,badTimes)
 
@@ -683,17 +683,80 @@ def QC(node):
                 elif ConfigFile.settings['SensorType'].lower() == 'trios':
                     Utilities.filterData(esGroup,badTimes,'L1AQC')
                     Utilities.filterData(liGroup,badTimes,'L1AQC')
-                    Utilities.filterData(ltGroup,badTimes,'L1AQC')                
+                    Utilities.filterData(ltGroup,badTimes,'L1AQC')
 
         # Next apply the Meteorological FLAGGING prior to slicing
-        esData = referenceGroup.getDataset("ES")        
+        # esData = referenceGroup.getDataset("ES")
         if enableMetQualityCheck:
             # msg = "Applying meteorological filtering to eliminate spectra."
             msg = "Applying meteorological flags. Met flags are NOT used to eliminate spectra."
             print(msg)
             Utilities.writeLogFile(msg)
             badTimes = ProcessL1bqc.metQualityCheck(referenceGroup, sasGroup, py6sGroup, ancGroup)
 
+        # NOTE: This is not finalized and needs a ConfigFile.setting #########################################
+        ConfigFile.settings['bL2filterMetFlags'] = 0
+        if ConfigFile.settings['bL2filterMetFlags'] == 1:
+            # Placeholder to filter out based on Met Flags
+            metFlags = ancGroup.datasets['MET_FLAGS']
+            AncDatetime = metFlags.columns['Datetime']
+            Flag3 = metFlags.columns['Flag3']
+
+            badTimes = []
+            for indx, dateTime in enumerate(AncDatetime):
+                if Flag3[indx]:
+                    badTimes.append(dateTime)
+
+            badTimes = np.unique(badTimes)
+            # Duplicate each element to a list of two elements in a list
+            badTimes = np.rot90(np.matlib.repmat(badTimes,2,1), 3)
+            # msg = f'{len(np.unique(badTimes))/len(AncDatetime)*100:.1f}% of spectra flagged'
+
+            if badTimes is not None:
+                msg = "Removing spectra from Met flags. ######################### Hard-coded override for Flag3"
+                print(msg)
+                Utilities.writeLogFile(msg)
+                check = Utilities.filterData(referenceGroup, badTimes)
+                if check > 0.99:
+                    msg = "Too few spectra remaining. Abort."
+                    print(msg)
+                    Utilities.writeLogFile(msg)
+                    return False
+                check = Utilities.filterData(sasGroup, badTimes)
+                if check > 0.99:
+                    msg = "Too few spectra remaining. Abort."
+                    print(msg)
+                    Utilities.writeLogFile(msg)
+                    return False
+                check = Utilities.filterData(ancGroup, badTimes)
+                if check > 0.99:
+                    msg = "Too few spectra remaining. Abort."
+                    print(msg)
+                    Utilities.writeLogFile(msg)
+                    return False
+
+                if py6sGroup is not None:
+                    Utilities.filterData(py6sGroup,badTimes)
+
+                # Filter L1AQC data for L1BQC criteria
+                if ConfigFile.settings['SensorType'].lower() == 'seabird':
+                    check = []
+                    check.append(Utilities.filterData(esDarkGroup,badTimes,'L1AQC'))
+                    check.append(Utilities.filterData(esLightGroup,badTimes,'L1AQC'))
+                    check.append(Utilities.filterData(liDarkGroup,badTimes,'L1AQC'))
+                    check.append(Utilities.filterData(liLightGroup,badTimes,'L1AQC'))
+                    check.append(Utilities.filterData(ltDarkGroup,badTimes,'L1AQC'))
+                    check.append(Utilities.filterData(ltLightGroup,badTimes,'L1AQC'))
+                    if any(np.array(check) > 0.99):
+                        msg = "Too few spectra remaining. Abort."
+                        print(msg)
+                        Utilities.writeLogFile(msg)
+                        return False
+                elif ConfigFile.settings['SensorType'].lower() == 'trios':
+                    Utilities.filterData(esGroup,badTimes,'L1AQC')
+                    Utilities.filterData(liGroup,badTimes,'L1AQC')
+                    Utilities.filterData(ltGroup,badTimes,'L1AQC')
+
         return True
 
     @staticmethod

Source/ProcessL2.py

@@ -2357,7 +2357,7 @@ def processL2(root,station=None):
                 newGrp = node.addGroup(grp.id)
                 newGrp.copy(grp)
                 for ds in newGrp.datasets:
-                    newGrp.datasets[ds].datasetToColumns()
+                    newGrp.datasets[ds].datasetToColumns()                    
 
         # Process stations, ensembles to reflectances, OC prods, etc.
         if not ProcessL2.stationsEnsemblesReflectance(node, root,station):