gnssrefl.sd_libs module

gnssrefl.sd_libs.RH_ortho_plot2(station, H0, year, txtdir, fs, time_rh, rh, gap_min_val, th, spline, delta_out, csvfile)

Makes a plot of the final spline fit to the data. Output time interval controlled by the user.

It also now writes out the file with the spline fit

  • station (str) – name of station, 4 ch

  • H0 (float) – datum correction (orthometric height) to convert RH to MSL data, in meters

  • year (int) – year of the time series (ultimately should not be needed)

  • txtdir (str) – location of plot

  • fs (int) – fontsize

  • time_rh (numpy of floats) – time of rh values, in fractional doy I believe

  • rh (numpy of floats) – refl hgt in meters

  • gap_min_val (float) – minimum length gap allowed, in day of year units

  • th (numpy floats) – time values in day of year units

  • spline (output of interpolate.BSpline) – used for fitting

  • delta_out (int) – how often spline is printed, in seconds

  • csvfile (bool) – print out csv instead of plain txt

gnssrefl.sd_libs.find_ortho_height(station, extension)

Find orthometric (sea level) height used in final subdaily spline output and plots. This value should be defined for the GPS L1 phase center of the GNSS antenna as this is what is assumed in the subdaily code.

  • station (str) – 4 ch station name

  • extension (str) – gnssir analysis, extension mode


Hortho – orthometric height from gnssir json analysis file as defined as Hortho, in meters. If your preferred value for Hortho is not present, it is calculated from the ellipsoidal height and EGM96.

Return type:


gnssrefl.sd_libs.mirror_plot(tnew, ynew, spl_x, spl_y, txtdir, station, beginT, endT)

Makes a plot of the spline fit to the mirrored RH data Plot is saved to txtdir as station_rhdot1.png

  • tnew (numpy of floats) – time in days of year, including the faked data, used for splines

  • ynew (numpy of floats) – RH in meters

  • spl_x (numpy of floats) – time in days of year

  • spl_y (numpy of floats) – smooth RH, meters

  • txtdir (str) – directory for plot

  • station (str) – name of station for title

  • beginT (float) – first time (day of year) real RH measurement

  • endT (float) – last time (day of year) for first real RH measurement


takes mjd array and converts to datetime for plotting.


mjd (numpy array of floats) – mod julian date



Return type:

numpy array of datetime objects

gnssrefl.sd_libs.numsats_plot(station, tval, nval, Gval, Rval, Eval, Cval, txtdir, fs, hires_figs, year)

makes the plot that summarizes the number of satellites in each constellation per epoch

  • station (str) – name of the station

  • tval (numpy array) – datetime objects?

  • nval (numpy array) – number of total satellites at a given epoch

  • Gval (numpy array) – number of galileo satellites at an epoch

  • Rval (numpy array) – number of glonass satellites at an epoch

  • Eval (numpy array) – number of galileo satellites at an epoch

  • Cval (numpy array) – number of beidou satellites at an epoch

  • txtdir (str) – where results are stored

  • fs (int) – fontsize for the plots

  • hires_figs (bool) – try to plot high resolution

  • year (int) – calendar year

gnssrefl.sd_libs.pickup_subdaily_json_defaults(xdir, station, extension)

picks up an existing gnssir analysis json. augments with subdaily parameters if needed. Returns the dictionary.

  • xdir (str) – REFL_CODE code location

  • station (str) – name of station

  • extension (str) – possible extension location


lsp – contents of gnssir json

Return type:


gnssrefl.sd_libs.print_badpoints(t, outliersize, txtdir, real_residuals)

prints outliers to a file.

  • t (numpy array) – lomb scargle result array of “bad points”. Format given below

  • outliersize (float) – outlier criterion, in meters

  • txtdir (str) – directory where file is written

  • real_residuals (numpy array of floats) – assume this is RH residuals in meters

Return type:

writes to a file called outliers.txt in the Files/station area

gnssrefl.sd_libs.quickTr(year, doy, frachours)

takes timing from lomb scargle code (year, doy) and UTC hour (fractional) and returns a date string

  • year (int) – full year

  • doy (int) – day of year

  • frachours (float) – real-valued UTC hour


datestring – date ala YYYY-MM-DD HH:MM:SS

Return type:


gnssrefl.sd_libs.rh_plots(otimes, tv, station, txtdir, year, d1, d2, percent99)

overview plots for rh_plot

  • otimes (numpy array of datetime objects) – observation times

  • tv (numpy array) – gnssrefl results written into this variable using loadtxt

  • station (str) – station name, only used for the title

  • txtdir (str) – directory where the plots will be written to

  • year (int) – what year is being analyzed

  • d1 (int) – minimum day of year

  • d2 (int) – maximum day of year

  • percent99 (bool) – whether you want only the 1-99 percentile plotted

gnssrefl.sd_libs.rhdot_plots(th, correction, rhdot_at_th, tvel, yvel, fs, station, txtdir, hires_figs, year)

makes the rhdot correction plots

  • th (numpy array) – time of obs, day of year

  • correction (numpy array) – rhcorrections in meters

  • rhdot_at_th (numpy array of floats) – spline fit for rhdot in meters

  • tvel (numpy array of floats) – time for surface velocity in days of year

  • yvel (numpy array of floats) – surface velocity in m/hr

  • fs (integer) – fontsize

  • station (str) – station name

  • txtdir (str) – file directory for output

  • hires_figs (bool) – whether you want eps instead of png

  • year (int) – calendar year

gnssrefl.sd_libs.stack_two_more(otimes, tv, ii, jj, stats, station, txtdir, sigma, kplt, hires_figs, year)

makes a plot of the reflector heights before and after minimal editing

  • otimes (numpy array of datetime objects) – observation times

  • tv (numpy array) – variable with the gnssrefl LSP results

  • ii (numpy array) – indices of good data

  • jj (numpy array) – indices of bad data

  • station (str) – station name

  • txtdir (str) – directory where plots will be written

  • sigma (float) – what kind of standard deviation is used for outliers (3sigma, 2.5 sigma etc)

  • kplt (bool) – make extra plot for kristine

  • year (int) – calendar year

gnssrefl.sd_libs.subdaily_resids_last_stage(station, year, th, biasCor_rh, spline_at_GPS, fs, strsig, hires_figs, txtdir, ii, jj, th_even, spline_whole_time)

Makes the final residual plot for subdaily (after RHdot and IF correction made). Returns the bad points …

Allows either the original or multiyear option..

  • station (str) – 4 character station name

  • year (int) – calendar year

  • th (numpy array of ??) – time variable of some kind, fractional day of year ?

  • biasCor_rh (numpy array of floats) – refl hgts that have been corrected for RHdot and IF

  • spline_at_GPS (numpy array of floats) – RH derived From the spline fit and calculated at GPS time tags

  • fs (int) – font size

  • strsig (str) – sigma string to go on the legend

  • hires_figs (bool) – whether to save the plots with better resolution

  • txtdir (str) – directory where the plot will be saved

  • ii (numpy array) – indices of the outliers?

  • jj (numpy array) – indices of the values to keep?

  • th_even (numpy array) – evenly spaced time values, day of year

  • spline_whole_time (numpy array of flots) – splinefit for ???


badpoints2 – RH residuals

Return type:

numpy array of floats

gnssrefl.sd_libs.testing_nvals(Gval, Rval, Eval, Cval)

writing the number of observations per constellation as a test. not currently used

Parameters Gval: numpy array

GPS RH values

Rvalnumpy array


Evalnumpy array

Galileo RH values


Beidou RH values

writes to a file - kristine.txt returns nothing

gnssrefl.sd_libs.two_stacked_plots(otimes, tv, station, txtdir, year, d1, d2, hires_figs)

This actually makes three stacked plots - not two, LOL It gives an overview for quality control

  • otimes (numpy array of datetime objects) – observations times

  • tv (numpy array) – gnssrefl results written into this variable using loadtxt

  • station (str) – station name, only used for the title

  • txtdir (str) – where the plots will be written to

  • year (int) – what year is being analyzed

  • d1 (int) – minimum day of year

  • d2 (int) – maximum day of year

  • hires_figs (bool) – true for eps instead of png

gnssrefl.sd_libs.write_spline_output(iyear, th, spline, delta_out, station, txtdir, Hortho)

Writing the output of the spline fit to the final RH time series. No output other than this text file for year 2023 and station name ssss:


I do not think this is used anymore. It has been consolidated with the plot code.

  • iyear (int) – full year

  • th (numpy array) – time values of some kind … maybe fractional day of years?

  • spline (fit, output of interpolate.BSpline) – needs doc

  • delta_out (int) – how often you want the splinefit water level written, in seconds

  • station (str) – station name

  • txtdir (str) – output directory

  • Hortho (float) – orthometric height used to convert RH to something more sea level like meters

gnssrefl.sd_libs.writejsonfile(ntv, station, outfile)

subdaily RH values written out in json format

This does not appear to be used

  • ntv (numpy of floats) – LSP results

  • station (str) – 4 ch station name

  • outfile (str) – filename for output

gnssrefl.sd_libs.writeout_spline_outliers(tvd_bad, txtdir, residual, filename)

Write splinefit outliers to a text file.

  • tvd_bad (numpy array) – output of the lomb scargle calculations

  • txtdir (str) – directory for the output, i.e. $REFL_CODE/FiLes/station

  • residual (numpy array) – RH outliers in units of meters (!)

  • filename (str) – name of file being written