allow some subdaily options to be stored in the gnssir_input created …

…json file

CHANGELOG.md

@@ -4,6 +4,12 @@ All notable changes to this project will be documented in this file.
 The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.0.0/).
 
 
+## 3.1.4 
+
+Added a gazillion options to gnssir_input to accommodate daily_avg and subdaily 
+users that want to store strategies there. See the gnssir_input code for more details,
+i.e. the parameter names.
+
 ## 3.1.3
 
 Allow users to store required input parameters to daily_avg in the json used

README.md

@@ -1,4 +1,4 @@
-# gnssrefl v3.1.3 
+# gnssrefl v3.1.4 
 
 [![PyPI Version](https://img.shields.io/pypi/v/gnssrefl.svg)](https://pypi.python.org/pypi/gnssrefl) 
 

developer_notes.md

@@ -120,3 +120,6 @@ Make sure people are using version 3.9
 is version 3.10 ok?  for sure 3.11 is not.
 
 How to make a new github release?
+
+myAnimation.gif is in the EarthScopeLecture folder in Presentations.
+Need to find a better place for it

gnssrefl/gnssir_input.py

@@ -36,20 +36,25 @@ def parse_arguments():
     parser.add_argument("-delTmax", default=None, type=float, help="max arc length (min) default is 75. Shorten for tides.")
     parser.add_argument("-frlist", nargs="*",type=int,  help="User defined frequencies using our nomenclature")
     parser.add_argument("-azlist2", nargs="*",type=float,  default=None,help="list of azimuth regions, default is 0-360") 
-    parser.add_argument("-ellist", nargs="*",type=float,  default=None,help="List of elevation angles to allow more complex analysis scenarios-advanced users only!") 
+    parser.add_argument("-ellist", nargs="*",type=float,  default=None,help="List of elevation angles, advanced users only!") 
     parser.add_argument("-refr_model", default=1, type=int, help="refraction model. default is 1, zero turns it off)")
     parser.add_argument("-Hortho", default=None, type=float, help="station orthometric height, meters")
     parser.add_argument("-pele", nargs="*", type=float, help="min and max elevation angle in direct signal removal, default is 5-30")
     parser.add_argument("-daily_avg_medfilter", default=None, type=float, help="daily_avg, median filter, meters")
     parser.add_argument("-daily_avg_reqtracks", default=None, type=int, help="daily_avg, ReqTracks parameter ")
-    parser.add_argument("-subdaily_knots", default=None, type=float, help="subdaily, knots")
-    parser.add_argument("-update", default=None, type=str, help="mode to update existing json (under dev)")
-
+    parser.add_argument("-subdaily_alt_sigma", default=None, type=str, help="subdaily, Nievinski sigma")
+    parser.add_argument("-subdaily_sigma", default=None, type=str, help="subdaily, sigma for part 1")
+    parser.add_argument("-subdaily_ampl", default=None, type=float, help="subdaily, LSP amplitude override")
+    parser.add_argument("-subdaily_delta_out", default=None, type=int, help="subdaily, spline output interval (sec)")
+    parser.add_argument("-subdaily_knots", default=None, type=float, help="subdaily, knots, how many per day")
+    parser.add_argument("-subdaily_subdir", default=None, type=str, help="subdaily, output directory")
+    parser.add_argument("-subdaily_spline_outlier1", default=None, type=float, help="subdaily, outlier value (m), part1")
+    parser.add_argument("-subdaily_spline_outlier2", default=None, type=float, help="subdaily, outlier value (m), part2")
 
     args = parser.parse_args().__dict__
 
     # convert all expected boolean inputs from strings to booleans
-    boolean_args = ['allfreq', 'l1', 'l2c', 'xyz', 'refraction','update']
+    boolean_args = ['allfreq', 'l1', 'l2c', 'xyz', 'refraction','subdaily_alt_sigma']
     args = str2bool(args, boolean_args)
 
     # only return a dictionary of arguments that were added from the user - all other defaults will be set in code below
@@ -61,7 +66,10 @@ def make_gnssir_input(station: str, lat: float=0, lon: float=0, height: float=0,
        ampl: float = 5.0, allfreq: bool = False, l1: bool = False, l2c: bool = False, 
        xyz: bool = False, refraction: bool = True, extension: str = '', ediff: float=2.0, 
        delTmax: float=75.0, frlist: list=[], azlist2: list=[0,360], ellist : list=[], refr_model : int=1, 
-                      Hortho : float = None, pele: list=[5,30], update: bool=False, daily_avg_reqtracks: int=None, daily_avg_medfilter=None):
+                      Hortho : float = None, pele: list=[5,30], daily_avg_reqtracks: int=None, 
+                      daily_avg_medfilter: float =None, subdaily_alt_sigma : bool=None, subdaily_ampl : float=None, subdaily_delta_out : float=None, 
+                      subdaily_knots : int=None, subdaily_sigma: float=None, subdaily_subdir: str=None, 
+                      subdaily_spline_outlier1: float=None, subdaily_spline_outlier2: float=None):
 
     """
     This new script sets the Lomb Scargle analysis strategy you will use in gnssir. It saves your inputs 
@@ -219,14 +227,36 @@ def make_gnssir_input(station: str, lat: float=0, lon: float=0, height: float=0,
     pele : float
         min and max elevation angles in direct signal removal, i.e. 3 40. Default is 5 30. 
 
-    update : bool
-        allows you to update an existing json
-
-    daily_avg_reqtracks : int
+    daily_avg_reqtracks : int, optional
         number of tracks required for daily_avg code
 
-    daily_avg_medfilter : float
+    daily_avg_medfilter : float, optional
         median filter value required for daily_avg code (meters)
+
+    subdaily_alt_sigma : bool, optional
+        use Nievinski sigma definition
+
+    subdaily_ampl : float, optional
+        override the required LSP amplitude  
+
+    subdaily_delta_out : int, optional
+        spacing for final subdaily spline output
+
+    subdaily_knots : int, optional
+        number of knots per day for subdaily spline fits
+
+    subdaily_sigma : float, optional
+        how many standard deviations for outliers in subdaily code setting
+
+    subdaily_subdir : str, optional
+        non-standard location for subdaily outputs
+
+    subdaily_spline_outlier1 : float, optional
+        alternate setting for outlier detection in part1
+
+    subdaily_spline_outlier2 : float, optional
+        alternate setting for outlier detection in part2
+
     """
 
     # make sure environment variables exist
@@ -237,10 +267,6 @@ def make_gnssir_input(station: str, lat: float=0, lon: float=0, height: float=0,
         print('station name must be four characters long. Exiting.')
         sys.exit()
 
-    if update:
-        print('You have selected to update an existing json')
-        lsp = guts2.read_json_file(station, extension)
-        sys.exit()
 
 # location of the site - does not have to be very good.  within 100 meters is fine
     query_unr = False
@@ -382,6 +408,7 @@ def make_gnssir_input(station: str, lat: float=0, lon: float=0, height: float=0,
     lsp['overwriteResults'] = True
 
     # if snr file does not exist, try to make one
+    # i don't think this option exists anymore
     lsp['seekRinex'] = False
 
     # compress snr files after analysis - saves disk space
@@ -404,18 +431,27 @@ def make_gnssir_input(station: str, lat: float=0, lon: float=0, height: float=0,
     lsp['delTmax'] = delTmax  
 
     # gzip SNR files after running the code
+    # this really should be set to True.  the code is obviously ignoring it
     lsp['gzip'] = False   
 
     lsp['ellist'] = ellist
 
     lsp['refr_model'] = refr_model
 
-    # added for people that want to save their daily average strategies.
+    # added for people that want to save their daily average and subdaily strategies.
     # if not set, then they are saved as None.
     lsp['daily_avg_reqtracks'] = daily_avg_reqtracks
-
     lsp['daily_avg_medfilter'] = daily_avg_medfilter
 
+    lsp['subdaily_alt_sigma'] = subdaily_alt_sigma
+    lsp['subdaily_ampl'] = subdaily_ampl
+    lsp['subdaily_delta_out'] = subdaily_delta_out
+    lsp['subdaily_knots'] = subdaily_knots
+    lsp['subdaily_sigma'] = subdaily_sigma
+    lsp['subdaily_spline_outlier1'] = subdaily_spline_outlier1
+    lsp['subdaily_spline_outlier2'] = subdaily_spline_outlier2
+    lsp['subdaily_subdir'] = subdaily_subdir
+
     print('writing out to:', outputfile)
     print(lsp)
     with open(outputfile, 'w+') as outfile:

gnssrefl/gps.py

@@ -6244,18 +6244,9 @@ def set_subdir(subdir):
         print('This will tell the code where to put the output.')
         sys.exit()
 
-    outdir = xdir  + '/Files/'
+    outdir = xdir  + '/Files/' + subdir 
     if not os.path.exists(outdir) :
-        subprocess.call(['mkdir', outdir])
-
-    if subdir == '':
-        okokk = 1
-        #print('Using this output directory: ', outdir)
-    else:
-        outdir = xdir  + '/Files/' + subdir + '/'
-        #print('Using this output directory: ', outdir)
-        if not os.path.exists(outdir) :
-            subprocess.call(['mkdir', outdir])
+        subprocess.call(['mkdir', '-p', outdir])
 
     return
 

gnssrefl/sd_libs.py

@@ -1041,6 +1041,7 @@ def RH_ortho_plot2( station, H0, year,  txtdir, fs, time_rh, rh, gap_min_val,th,
 
     # find indices of gaps  that are larger than a certain value
     ii = (gdiff > gap_min_val)
+    #print('gap_min_val',gap_min_val, gdiff)
     N = len(mjd[ii])
     Ngdiff = len(gdiff)
 
@@ -1087,3 +1088,58 @@ def RH_ortho_plot2( station, H0, year,  txtdir, fs, time_rh, rh, gap_min_val,th,
     g.save_plot(pfile)
 
     return
+
+def pickup_subdaily_json_defaults(xdir, station, extension):
+    """
+    picks up an existing gnssir analysis json. augments
+    with subdaily parameters if needed. Returns the dictionary.
+
+    Parameters
+    ----------
+    xdir : str
+        REFL_CODE code location
+    station : str
+        name of station
+    extension : str
+        possible extension location
+
+    Returns
+    -------
+    lsp : dictionary 
+        contents of gnssir json
+    """
+
+    # changed default in subdaily to None instead of '' So have to check
+    if extension is None:
+        lsp = guts2.read_json_file(station, '',noexit=True)
+    else:
+        lsp = guts2.read_json_file(station, extension,noexit=True)
+
+    # check json for subdaily settings
+    # if values not found, then set them to None
+
+    if 'subdaily_ampl' not in lsp:
+        lsp['subdaily_ampl'] = None
+
+    if 'subdaily_spline_outlier1' not in lsp:
+        lsp['subdaily_spline_outlier1'] = None
+
+    if 'subdaily_spline_outlier2' not in lsp:
+        lsp['subdaily_spline_outlier2'] = None
+
+    if 'subdaily_knots' not in lsp:
+        lsp['subdaily_knots'] = None
+
+    if 'subdaily_delta_out' not in lsp:
+        lsp['subdaily_delta_out'] = None
+
+    if 'subdaily_alt_sigma' not in lsp:
+        lsp['subdaily_alt_sigma'] = None
+
+    if 'subdaily_sigma' not in lsp:
+        lsp['subdaily_sigma'] = None
+
+    if 'subdaily_subdir' not in lsp:
+        lsp['subdaily_subdir'] = None
+
+    return lsp

gnssrefl/subdaily.py

@@ -786,6 +786,7 @@ def rhdot_correction2(station,fname,fname_new,pltit,outlierV,outlierV2,**kwargs)
     # output will go to REFL_CODE/Files unless txtdir provided
     xdir = os.environ['REFL_CODE']
 
+    # change units to day of year
     gap_min_val = kwargs.get('gap_min_val',6.0)
     gap_min_val = gap_min_val/24 # change to DOY units
 
@@ -1128,7 +1129,7 @@ def rhdot_correction2(station,fname,fname_new,pltit,outlierV,outlierV2,**kwargs)
         tnew, ynew = flipit(tvd_new,column)
 
     Ndays = tnew.max()-tnew.min()
-    print('trying knots2')
+    #print('trying knots2')
     knots_per_day = knots2_per_day
 
     numKnots = int(knots_per_day*(Ndays))
@@ -1188,7 +1189,6 @@ def rhdot_correction2(station,fname,fname_new,pltit,outlierV,outlierV2,**kwargs)
     # write out the files with RH dot and IF bias corrected - and again without 3 sigma outliers
     bias_corrected_filename = fname_new + 'IF'; extraline = ''; 
     biasCor_rh = tvd_new[jj,24]
-    print(writecsv)
     write_subdaily(bias_corrected_filename,station,tvd_new[jj,:], writecsv,extraline, newRH_IF=biasCor_rh)
 
     new_outliers = tvd_new[ii,:]