diff --git a/.idea/inspectionProfiles/profiles_settings.xml b/.idea/inspectionProfiles/profiles_settings.xml new file mode 100644 index 0000000..105ce2d --- /dev/null +++ b/.idea/inspectionProfiles/profiles_settings.xml @@ -0,0 +1,6 @@ + + + + \ No newline at end of file diff --git a/.idea/misc.xml b/.idea/misc.xml deleted file mode 100644 index 789e113..0000000 --- a/.idea/misc.xml +++ /dev/null @@ -1,8 +0,0 @@ - - - - - - - diff --git a/RunWorkFlow_STWAVE.py b/RunWorkFlow_STWAVE.py index 84cafa7..f0c0c15 100755 --- a/RunWorkFlow_STWAVE.py +++ b/RunWorkFlow_STWAVE.py @@ -8,6 +8,7 @@ import os, getopt, sys, shutil, glob, platform, logging, yaml from testbedutils import fileHandling from getdatatestbed import getDataFRF +import pickle def Master_STWAVE_run(inputDict): """This will run STWAVE with any version prefix given start, end, and timestep @@ -81,7 +82,7 @@ def Master_STWAVE_run(inputDict): # ______________________________gather all data _____________________________ if generateFlag == True: go = getDataFRF.getObs(projectStart, projectEnd) # initialize get observation - rawspec = go.getWaveSpec(gaugenumber='waverider-26m', specOnly=True) + rawspec = go.getWaveData(gaugenumber='waverider-26m', spec=True) rawWL = go.getWL() rawwind = go.getWind(gaugenumber=0) loc_dict = go.get_sensor_locations(datafile=FRFgaugelocsFile, window_days=14) @@ -92,18 +93,20 @@ def Master_STWAVE_run(inputDict): # run the process through each of the above dates errors, errorDates, curdir = [], [], os.getcwd() for time in dateStringList: - print(' ------------------------------ START %s --------------------------------' %time) - + # print(' ------------------------------ START %s --------------------------------' %time) + # try: datadir = os.path.join(workingDirectory, ''.join(time.split(':'))) # moving to the new simulation's folder - if generateFlag == True: - [nproc_par, nproc_nest] = STsimSetup(time, inputDict, rawwind, rawWL, rawspec, bathy, loc_dict) + pickleSaveFname = os.path.join(datadir, ''.join(time.split(':')) + '_io.pickle') - if nproc_par == -1 or nproc_nest == -1: - print('************************\nNo Data available\naborting run\n***********************') - # remove generated files? - shutil.rmtree(os.path.join(workingDirectory,''.join(time.split(':')))) - continue # this is to return to the next time step if there's no data + if generateFlag == True: + nproc_par, nproc_nest, stio = STsimSetup(time, inputDict, rawwind, rawWL, rawspec, bathy) + pickle.dump(stio, open(pickleSaveFname, 'wb'), protocol=pickle.HIGHEST_PROTOCOL) + if nproc_par == -1 or nproc_nest == -1: + print('************************\nNo Data available\naborting run\n***********************') + # remove generated files? + shutil.rmtree(os.path.join(workingDirectory,''.join(time.split(':')))) + continue # this is to return to the next time step if there's no data if runFlag == True: os.chdir(datadir) # changing locations to where simulation files live @@ -127,10 +130,13 @@ def Master_STWAVE_run(inputDict): count = nproc_nest # lower the processors called for to match sim file (otherwise will throw segfault) child = check_output('mpiexec -n {} {} {}nested.sim'.format(count, executableLocation, ''.join(time.split(':'))), shell=True) print((' Simulations took {:.2f} hours'.format((DT.datetime.now() - t).total_seconds()/3600))) + pickle.dump(stio, open(pickleSaveFname, 'wb'), protocol=pickle.HIGHEST_PROTOCOL) # run analyze and archive script os.chdir(curdir) # change back after runing simulation locally if analyzeFlag == True: - beachWaves = STanalyze(time, inputDict) + if runFlag is False: + stio = pickle.load(open(pickleSaveFname, 'rb')) + beachWaves = STanalyze(time, inputDict, stio) if pFlag == True and DT.date.today() == endTime.date(): print('**\n Moving Plots! \n &&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&') # move files diff --git a/frontback/frontBackFUNWAVE.py b/frontback/frontBackFUNWAVE.py index 0d3aa72..772a836 100644 --- a/frontback/frontBackFUNWAVE.py +++ b/frontback/frontBackFUNWAVE.py @@ -16,6 +16,7 @@ from testbedutils import sblib as sb from testbedutils import waveLib as sbwave from testbedutils import fileHandling +from testbedutils import timeDomain from plotting.operationalPlots import obs_V_mod_TS from testbedutils import geoprocess as gp import multiprocessing @@ -63,14 +64,6 @@ def FunwaveSimSetup(startTime, rawWL, rawspec, bathy, inputDict): assert 'time' in rawspec, "\n++++\nThere's No Wave data" # preprocess wave spectra - - #if version_prefix.lower() == 'base': - # wavepacket1 = prepdata.prep_SWASH_spec(rawspec, version_prefix, model=model, nf=inputDict['modelSettings']['nf']) - - #else: - # #raise NotImplementedError('pre-process TS data ') - # wavepacket1 = prepdata.prep_SWASH_spec(rawspec, version_prefix, model=model, nf=inputDict['modelSettings']['nf']) - wavepacket = prepdata.prep_SWASH_spec(rawspec, version_prefix, model=model, nf=nf, phases=phases, grid=inputDict['modelSettings']['grid']) @@ -185,8 +178,8 @@ def FunwaveAnalyze(startTime, inputDict, fio): ###################################################################################################################### ###################################################################################################################### - outputFolder = os.path.join(fpath,'output') - print('Loading files ',outputFolder) + outputFolder = os.path.join(fpath, 'output') + print('Loading files ', outputFolder) ## upload depth file depthFile = os.path.join(outputFolder, 'dep.out') @@ -216,34 +209,36 @@ def FunwaveAnalyze(startTime, inputDict, fio): eta = simData['eta'].squeeze() # now adapting Chuan's runup code, here we use 0.08 m for runup threshold - r_depth = 0.08 # 4.0 * np.nanmax(np.abs(h[runupInd][1:] - h[runupInd][:-1])) - + # r_depth = 0.08 # 4.0 * np.nanmax(np.abs(h[runupInd][1:] - h[runupInd][:-1])) # Preallocate runup variable - runup = np.zeros(eta.shape[0]) - x_runup = np.zeros_like(runup) - - for aa in range(runup.shape[0]): - # Water depth - wdepth = eta[aa, :] + simData['elevation'] - # Find the runup contour (search from left to right) - wdepth_ind = np.argmin(abs(wdepth - r_depth)) # changed from Chuan's original code - # Store the water surface elevation in matrix - runup[aa] = eta[aa, wdepth_ind] # unrealistic values for large r_depth - # runup[aa]= -h[wdepth_ind] - # Store runup position - x_runup[aa] = simData['xFRF'][wdepth_ind] - maxRunup = np.amax(runup) - + # runup = np.zeros(eta.shape[0]) + # x_runup = np.zeros_like(runup) + # + # for aa in range(runup.shape[0]): + # # Water depth + # wdepth = eta[aa, :] + simData['elevation'] + # # Find the runup contour (search from left to right) + # wdepth_ind = np.argmin(abs(wdepth - r_depth)) # changed from Chuan's original code + # # Store the water surface elevation in matrix + # runup[aa] = eta[aa, wdepth_ind] # unrealistic values for large r_depth + # # runup[aa]= -h[wdepth_ind] + # # Store runup position + # x_runup[aa] = simData['xFRF'][wdepth_ind] + # maxRunup = np.amax(runup) + + runupTS, x_runup, r_depth = timeDomain.runup_func(eta, Depth1D, simData['xFRF'], r_depth=0.1) + r2, peaks, maxSetup = timeDomain.identifyR2(runupTS, percentile=2) + r2 = r2 + simMeta['WL'] ###################################################################################################################### ###################################################################################################################### ################################## plotting ######################################################################### ###################################################################################################################### ###################################################################################################################### - fileHandling.makeCMTBfileStructure(path_prefix,date_str=datestring) + fileHandling.makeCMTBfileStructure(path_prefix, date_str=datestring) figureBaseFname = 'CMTB_waveModels_{}_{}_'.format(model, version_prefix) # make function for processing timeseries data - data = simData['eta'].squeeze()[cutRampingTime:,:] + data = simData['eta'].squeeze()[cutRampingTime:, :] time = [] for i in range(len(simData['time'].squeeze()[cutRampingTime:])): ## change time from float to datetime @@ -253,7 +248,7 @@ def FunwaveAnalyze(startTime, inputDict, fio): SeaSwellCutoff = 0.05 # cutoff between sea/swell and IG nSubSample = 5 - fspec, freqs = sbwave.timeSeriesAnalysis1D(np.asarray(time), data, bandAvg=3)#6,WindowLength=20) + fspec, freqs = sbwave.timeSeriesAnalysis1D(np.asarray(time), data, bandAvg=3) #6,WindowLength=20) total = sbwave.stats1D(fspec=fspec, frqbins=freqs, lowFreq=None, highFreq=None) SeaSwellStats = sbwave.stats1D(fspec=fspec, frqbins=freqs, lowFreq=SeaSwellCutoff, highFreq=None) IGstats = sbwave.stats1D(fspec=fspec, frqbins=freqs, lowFreq=None, highFreq=SeaSwellCutoff) @@ -264,6 +259,7 @@ def FunwaveAnalyze(startTime, inputDict, fio): ############################################################################################################# WL = simMeta['WL'] #added in editing, should possibly be changed? setup = np.mean(simData['eta'] + WL, axis=0).squeeze() + if plotFlag == True: from plotting import operationalPlots as oP ## remove images before making them if reprocessing @@ -329,8 +325,8 @@ def FunwaveAnalyze(startTime, inputDict, fio): 'waveHsIG': np.expand_dims(IGstats['Hm0'], axis=0), 'elevation': np.expand_dims(simData['elevation'], axis=0), 'eta': np.expand_dims(simData['eta'], axis=0), - 'totalWaterLevel': maxRunup, - 'totalWaterLevelTS': np.expand_dims(runup, axis=0), + 'totalWaterLevel': r2, + 'totalWaterLevelTS': np.expand_dims(runupTS, axis=0), 'velocityU': np.expand_dims(simData['velocityU'], axis=0), 'velocityV': np.expand_dims(simData['velocityV'], axis=0), 'waveHs': np.expand_dims(SeaSwellStats['Hm0'], axis=0), # or from HsTS?? diff --git a/frontback/frontBackSTWAVE.py b/frontback/frontBackSTWAVE.py index c281853..7fdd657 100755 --- a/frontback/frontBackSTWAVE.py +++ b/frontback/frontBackSTWAVE.py @@ -22,7 +22,7 @@ from testbedutils import fileHandling from subprocess import check_output -def STsimSetup(startTime, inputDict,allWind , allWL, allWave, bathy, loc_dict=None): +def STsimSetup(startTime, inputDict,allWind , allWL, allWave, bathy): """This Function is the master call for the data preparation for the Coastal Model STWAVE runs. It is designed to be handed data then utilize prep_datalib for model pre-processing. All Files are labeled by @@ -34,7 +34,6 @@ def STsimSetup(startTime, inputDict,allWind , allWL, allWave, bathy, loc_dict=No allWind(dict): from get data with all wind from entire project set allWL(dict): from get data with all waterlevel from entire project set allWave(dict): from get data with all wave from entire project set - gaugelocs(list): provides save points (lat/lon) Returns: nproc_parent (int): number of processors to run simultion for paret sim, will return -1 if this @@ -92,7 +91,8 @@ def STsimSetup(startTime, inputDict,allWind , allWL, allWave, bathy, loc_dict=No # __________________initalize needed classes ___________________________________ prepdata = STPD.PrepDataTools() - stio = inputOutput.stwaveIO('') # initializing io here so grid text can be written out + stio = inputOutput.stwaveIO(os.path.join(path_prefix, dateString)) # initializing io here so grid + # text can be written out ################################################################################################################### ####################### Begin Gathering Data ############################################################### @@ -105,7 +105,7 @@ def STsimSetup(startTime, inputDict,allWind , allWL, allWave, bathy, loc_dict=No # rotate, time match, interpolate wave spectra waveTimeList = sb.createDateList(d1, d2-DT.timedelta(seconds=3600), DT.timedelta(seconds=3600)) # time match list wavepacket = prepdata.prep_spec(rawspec, version_prefix, datestr=dateString, plot=plotFlag, full=full, - outputPath=path_prefix, waveTimeList=waveTimeList) + outputPath=os.path.join(path_prefix, dateString), waveTimeList=waveTimeList) print("number of wave records %d with %d interpolated points" % (np.shape(wavepacket['spec2d'])[0], sum(wavepacket['flag']))) # ____________ BATHY ______________________ @@ -158,25 +158,18 @@ def STsimSetup(startTime, inputDict,allWind , allWL, allWave, bathy, loc_dict=No print(' TODO: handle TDS location grabber') dataLocations = query(d1, d2, inputName='/home/spike/repos/TDSlocationGrabber/database', type='waves') # # get gauge nodes x/y new idea: put gauges into input/output instance for the model, then we can save it - statloc = [] + savePoints, statloc = [], [] for _, gauge in enumerate(['waverider-26m', 'waverider-17m', 'awac-11m', '8m-array', 'awac-6m', 'awac-4.5m', 'adop-3.5m', 'xp200m', 'xp150m', 'xp125m']): ii = np.argwhere(gauge == dataLocations['Sensor']).squeeze() - coord = gp.FRFcoord(dataLocations['Lon'][ii], dataLocations['Lat'][ii], coordType='LL') - print(' save point added: {} at xFRF {:.1f} yFRF {:.1f}'.format(gauge, coord['xFRF'], coord['yFRF'])) - statloc.append([coord['StateplaneE'], coord['StateplaneN']]) + if np.size(ii) > 0: + coord = gp.FRFcoord(dataLocations['Lon'][ii], dataLocations['Lat'][ii], coordType='LL') + print(' save point added: {} at xFRF {:.1f} yFRF {:.1f}'.format(gauge, coord['xFRF'], coord['yFRF'])) + statloc.append([coord['StateplaneE'], coord['StateplaneN']]) + savePoints.append(gauge) statloc = np.array(statloc) - - #old way to grab gauge data - # statloc = [] - # for gauge in list(loc_dict.keys()): - # coords = loc_dict[gauge] - # try: - # statloc.append([coords['spE'], coords['spN']]) - # except KeyError: - # continue - # statloc = np.array(statloc) - + stio.savePoints = savePoints + # assign nesting points in i/j grid coordinates if runNested is True: @@ -184,7 +177,7 @@ def STsimSetup(startTime, inputDict,allWind , allWL, allWave, bathy, loc_dict=No whichpointsy = np.linspace(gridNodesNested['yFRF'][0], gridNodesNested['yFRF'][-1], numNest) nestLocDict = {} #initalize nesting output locations for key in whichpointsy: - nestLocDict[key]= {'xFRF': gridNodesNested['xFRF'].max(),'yFRF':key} + nestLocDict[key] = {'xFRF': gridNodesNested['xFRF'].max(),'yFRF':key} else: nestLocDict = False ################################################################################################################### @@ -206,9 +199,10 @@ def STsimSetup(startTime, inputDict,allWind , allWL, allWave, bathy, loc_dict=No statloc=statloc, full=full) else: nproc_nest = None - return nproc_parent, nproc_nest + + return nproc_parent, nproc_nest, stio -def STanalyze(startTime, inputDict): +def STanalyze(startTime, inputDict, stio): """the master call for the simulation post process for STWAVE plots and netCDF files are made at request @@ -221,6 +215,8 @@ def STanalyze(startTime, inputDict): 'key' startTime (str): a string that has date in it by which the end of the run is designated ex: '2015-12-25T00:00:00Z' + stio(obcect): initalized class for model setup (with metadata and save points) + Returns: None @@ -228,7 +224,7 @@ def STanalyze(startTime, inputDict): # ___________________ Unpack input dictionary _________________________________ plotFlag = inputDict.get('plotFlag', True) model = inputDict.get('model', 'stwave') - version_prefix = inputDict['modelSettings']['version_prefix'] + version_prefix = inputDict['modelSettings']['version_prefix'].lower() path_prefix = inputDict.get('path_prefix', "{}".format(version_prefix)) simulationDuration = inputDict['simulationDuration'] Thredds_Base = inputDict.get('netCDFdir', '/home/{}/thredds_data/'.format(check_output('whoami', shell=True)[:-1])) @@ -256,7 +252,7 @@ def STanalyze(startTime, inputDict): ################################################################################################# ################################################################################################# - stio = stwaveIO(fpath) # =pathbase) looks for model output files in folder to analyze + # stios = stioR.stwaveIO(fpath) # =pathbase) looks for model output files in folder to analyze assert len(stio.obsefname) != 0, 'There are no data to load for %s' %startTime print('Loading Statistic Files ....') d = DT.datetime.now() # for data load timing @@ -279,7 +275,7 @@ def STanalyze(startTime, inputDict): # correct angles modelpacket_nest['WaveDm'] = anglesLib.angle_correct(modelpacket_nest['WaveDm']) modelpacket_nest['Udir'] = anglesLib.angle_correct(modelpacket_nest['Udir']) # wind direction - except IndexError: + except (IndexError, FileNotFoundError): nest = 0 # Load Spectral Data sets obse_packet = stio.obseload(nested=False) @@ -319,14 +315,12 @@ def STanalyze(startTime, inputDict): print(' ..begin loading spatial files ....') dep_pack = prepdata.GetOriginalGridFromSTWAVE(stio.simfname[0], stio.depfname[0]) Tp_pack = stio.genLoad(nested=False, key='waveTp') - # Tp_pack = stio.TPload(nested=0) # this function is currently faster than genLoad (a generalized load function) wave_pack = stio.genLoad(nested=False, key='wave') if runNested: rad_nest = stio.genLoad('rad', nested=True) break_nest = stio.genLoad('break', nested=True) dep_nest = prepdata.GetOriginalGridFromSTWAVE(stio.simfname_nest[0], stio.depfname_nest[0]) Tp_nest = stio.genLoad(key='waveTp', nested=True) - # wave_pack2 = stio.genLoad('wave', nested=False) wave_nest = stio.genLoad(key='wave', nested=True) print('Files loaded, in %s ' %(DT.datetime.now() - d)) @@ -511,8 +505,8 @@ def STanalyze(startTime, inputDict): cmtbLocalFldrArch = os.path.join(model, version_prefix) varYml = 'yaml_files/waveModels/{}/Field_var.yml'.format(model) - regGlobYml = 'yaml_files/waveModels/{}/{}/Field_Regional_{}_globalmeta.yml'.format(model, version_prefix, - version_prefix) + regGlobYml = 'yaml_files/waveModels/{}/{}/Field_Regional_{}_globalmeta.yml'.format(model, version_prefix.lower(), + version_prefix.lower()) outFileName = fileHandling.makeTDSfileStructure(Thredds_Base, cmtbLocalFldrArch, dateString, 'Regional-Field') assert os.path.isfile(regGlobYml), 'NetCDF yaml files are not created' makenc.makenc_field(data_lib=regionalDataLib, globalyaml_fname=regGlobYml, flagfname=flagfname, @@ -573,7 +567,7 @@ def STanalyze(startTime, inputDict): RegionalStations = ['waverider-26m', 'waverider-17m', 'awac-11m', '8m-array', 'awac-6m', 'awac-4.5m', 'adop-3.5m', 'xp200m', 'xp150m', 'xp125m'] # where :'s are for sims in the nested # writing station files from regional/parent simulation - for gg, station in enumerate(RegionalStations): + for gg, station in enumerate(stio.savePoints): if station != ':': # getting lat lon coords = gp.FRFcoord(stat_packet['Easting'][0, gg], stat_packet['Northing'][0, gg]) @@ -585,13 +579,13 @@ def STanalyze(startTime, inputDict): 'waveTp': stat_packet['Tp'][:, gg], 'waterLevel': stat_packet['WL'][:, gg], 'Umag': stat_packet['Umag'][:, gg], - 'Udir': stat_packet['Udir'][:, gg ], - 'Northing' : stat_packet['Northing'][0, gg], + 'Udir': stat_packet['Udir'][:, gg], + 'Northing': stat_packet['Northing'][0, gg], 'Easting': stat_packet['Easting'][0, gg], - 'Latitude' : coords['Lat'], - 'Longitude' : coords['Lon'], + 'Latitude': coords['Lat'], + 'Longitude': coords['Lon'], 'station_name': station, - 'directionalWaveEnergyDensity': obse_packet['ncSpec'][:,gg,:,:], + 'directionalWaveEnergyDensity': obse_packet['ncSpec'][:, gg, :, :], 'waveDirectionBins': obse_packet['ncDirs'], 'waveFrequency': obse_packet['Frequencies'], 'DX': dep_pack['DX'], @@ -613,90 +607,92 @@ def STanalyze(startTime, inputDict): # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # ## COLLECT ALL data (process, rotate, time pair and make comparison plots) # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # - go=getDataFRF.getObs(startTime, endTime) - if plotFlag == True and np.size(stat_packet['time']) > 1: - from testbedutils import waveLib as sbwave - print(' Plotting Time Series Data ') - stationList = ['waverider-26m', 'waverider-17m', 'awac-11m', '8m-array', 'awac-6m', 'awac-4.5m', 'adop-3.5m', 'xp200m', 'xp150m', 'xp125m'] - for gg, station in enumerate(stationList): - print('working on %s' %station) - # go get comparison data - w = go.getWaveSpec(station) - if w is not None and 'time' in w: # if there's data (not only location) - if station in go.directionalWaveGaugeList: - if full == False and station in go.directionalWaveGaugeList: - w['dWED'], w['wavedirbin'] = sbwave.HPchop_spec(w['dWED'], w['wavedirbin'], angadj=angadj) - obsStats = sbwave.waveStat(w['dWED'], w['wavefreqbin'], w['wavedirbin']) - else: # calc non directionalWaveGaugeList stats - obsStats = sbwave.stats1D(w['fspec'], w['wavefreqbin']) - - if station in ['waverider-17m', 'awac-11m', 'waverider-26m']: - modStats = sbwave.waveStat(obse_packet['ncSpec'][:, gg, :, :], obse_packet['Frequencies'], - obse_packet['ncDirs']) # compute model stats here - else: - if runNested is True: - modStats = sbwave.waveStat(obse_nested['ncSpec'][:, gg, :, :], obse_nested['Frequencies'], - obse_nested['ncDirs']) # compute model stats here - - # time match data - time, obsi, modi = sb.timeMatch(w['epochtime'], - np.arange(w['time'].shape[0]), - nc.date2num(stat_packet['time'][:], 'seconds since 1970-01-01'), - np.arange(len(stat_packet['time']))) # time match - # don't plot if theres only 1 dot on the plot... save time - if station in go.directionalWaveGaugeList: - plotList = ['Hm0', 'Tm', 'sprdF', 'sprdD', 'Tp', 'Dm'] - else: - plotList = ['Hm0', 'Tm', 'sprdF', 'Tp'] - for param in modStats: - if np.size(obsi) > 1 and np.size(modi) > 1 and param in plotList: - if param in ['Tm', 'Tp', 'Tm10', 'Tave']: - units = 's' - title = '%s period' % param - elif param in ['Hm0']: - units = 'm' - title = 'Wave Height %s ' %param - elif param in ['Dm', 'Dmp', 'Dp']: - units = 'degrees' - title = 'Direction %s' %param - elif param in ['sprdF']: - title = 'Spread %s' % param - units = 'Hz' - elif param in ['sprdD']: - units = 'degrees' - title = 'Spread %s ' % param - - ofname = os.path.join(path_prefix, dateString, - 'figures/CMTB-waveModels_{}_{}_station_{}_{}_{}.png'.format(model, version_prefix, - station, param, dateString)) - print('plotting ' + ofname) - dataDict = {'time': nc.num2date(time, 'seconds since 1970-01-01', - only_use_cftime_datetimes=False), - 'obs': obsStats[param][obsi.astype(np.int)], - 'model': modStats[param][modi.astype(np.int)], - 'var_name': param, - 'units': units, - 'p_title': title} - oP.obs_V_mod_TS(ofname, dataDict, logo_path='ArchiveFolder/CHL_logo.png') - if station == 'waverider-26m' and param == 'Hm0': - print(' skipping boundary spectral comparison') - continue - # this is a fail safe to abort run if the boundary conditions don't - # meet quality standards below - bias = 0.1 # bias has to be within 10 centimeters - RMSE = 0.1 # RMSE has to be within 10 centimeters - if isinstance(dataDict['obs'], np.ma.masked_array) and ~dataDict['obs'].mask.any(): - dataDict['obs'] = np.array(dataDict['obs']) - stats = sb.statsBryant(dataDict['obs'], dataDict['model']) - try: - assert stats['RMSE'] < RMSE, 'RMSE test on spectral boundary energy failed' - assert np.abs(stats['bias']) < bias, 'bias test on spectral boundary energy failed' - except: - print('!!!!!!!!!!FAILED BOUNDARY!!!!!!!!') - print('deleting data from thredds!') - # os.remove(fieldOfname) - os.remove(outFileName) - raise RuntimeError('The Model Is not validating its offshore boundary condition') + # go = getDataFRF.getObs(startTime, endTime) + # if plotFlag == True and np.size(stat_packet['time']) > 1: + # from testbedutils import waveLib as sbwave + # print(' Plotting Time Series Data ') + # stationList = stio.savePoints + # #['waverider-26m', 'waverider-17m', 'awac-11m', '8m-array', 'awac-6m', 'awac-4.5m', + # # 'adop-3.5m', 'xp200m', 'xp150m', 'xp125m'] + # for gg, station in enumerate(stationList): + # print('working on %s' %station) + # # go get comparison data + # w = go.getWaveData(station, spec=True) + # if w is not None and 'time' in w: # if there's data (not only location) + # if station in go.directionalWaveGaugeList: + # if full == False and station in go.directionalWaveGaugeList: + # w['dWED'], w['wavedirbin'] = sbwave.HPchop_spec(w['dWED'], w['wavedirbin'], angadj=angadj) + # obsStats = sbwave.waveStat(w['dWED'], w['wavefreqbin'], w['wavedirbin']) + # else: # calc non directionalWaveGaugeList stats + # obsStats = sbwave.stats1D(w['fspec'], w['wavefreqbin']) + # + # if station in ['waverider-17m', 'awac-11m', 'waverider-26m']: + # modStats = sbwave.waveStat(obse_packet['ncSpec'][:, gg, :, :], obse_packet['Frequencies'], + # obse_packet['ncDirs']) # compute model stats here + # else: + # if runNested is True: + # modStats = sbwave.waveStat(obse_nested['ncSpec'][:, gg, :, :], obse_nested['Frequencies'], + # obse_nested['ncDirs']) # compute model stats here + # + # # time match data + # time, obsi, modi = sb.timeMatch(w['epochtime'], + # np.arange(w['time'].shape[0]), + # nc.date2num(stat_packet['time'][:], 'seconds since 1970-01-01'), + # np.arange(len(stat_packet['time']))) # time match + # # don't plot if theres only 1 dot on the plot... save time + # if station in go.directionalWaveGaugeList: + # plotList = ['Hm0', 'Tm', 'sprdF', 'sprdD', 'Tp', 'Dm'] + # else: + # plotList = ['Hm0', 'Tm', 'sprdF', 'Tp'] + # for param in modStats: + # if np.size(obsi) > 1 and np.size(modi) > 1 and param in plotList: + # if param in ['Tm', 'Tp', 'Tm10', 'Tave']: + # units = 's' + # title = '%s period' % param + # elif param in ['Hm0']: + # units = 'm' + # title = 'Wave Height %s ' %param + # elif param in ['Dm', 'Dmp', 'Dp']: + # units = 'degrees' + # title = 'Direction %s' %param + # elif param in ['sprdF']: + # title = 'Spread %s' % param + # units = 'Hz' + # elif param in ['sprdD']: + # units = 'degrees' + # title = 'Spread %s ' % param + # + # ofname = os.path.join(path_prefix, dateString, + # 'figures/CMTB-waveModels_{}_{}_station_{}_{}_{}.png'.format(model, version_prefix, + # station, param, dateString)) + # print('plotting ' + ofname) + # dataDict = {'time': nc.num2date(time, 'seconds since 1970-01-01', + # only_use_cftime_datetimes=False), + # 'obs': obsStats[param][obsi.astype(np.int)], + # 'model': modStats[param][modi.astype(np.int)], + # 'var_name': param, + # 'units': units, + # 'p_title': title} + # oP.obs_V_mod_TS(ofname, dataDict, logo_path='ArchiveFolder/CHL_logo.png') + # if station == 'waverider-26m' and param == 'Hm0': + # print(' skipping boundary spectral comparison') + # continue + # # this is a fail safe to abort run if the boundary conditions don't + # # meet quality standards below + # bias = 0.1 # bias has to be within 10 centimeters + # RMSE = 0.1 # RMSE has to be within 10 centimeters + # if isinstance(dataDict['obs'], np.ma.masked_array) and ~dataDict['obs'].mask.any(): + # dataDict['obs'] = np.array(dataDict['obs']) + # stats = sb.statsBryant(dataDict['obs'], dataDict['model']) + # try: + # assert stats['RMSE'] < RMSE, 'RMSE test on spectral boundary energy failed' + # assert np.abs(stats['bias']) < bias, 'bias test on spectral boundary energy failed' + # except: + # print('!!!!!!!!!!FAILED BOUNDARY!!!!!!!!') + # print('deleting data from thredds!') + # # os.remove(fieldOfname) + # os.remove(outFileName) + # raise RuntimeError('The Model Is not validating its offshore boundary condition') # writing data out for UQ effort. Starting with parent domain (can complicate further as necessary) outData = regionalDataLib return outData diff --git a/frontback/frontBackSWASH.py b/frontback/frontBackSWASH.py index 7cf37d8..26fcee7 100644 --- a/frontback/frontBackSWASH.py +++ b/frontback/frontBackSWASH.py @@ -95,7 +95,7 @@ def SwashSimSetup(startTime, inputDict): except (RuntimeError, TypeError): WLpacket = None ### ____________ Get bathy grid from thredds ________________ - bathy = gdTB.getBathyIntegratedTransect(method=1, ybound=[940, 950]) + bathy = gdTB.getBathyIntegratedTransect(method=1, ybounds=[940, 950]) swsinfo, gridDict = prepdata.prep_SwashBathy(wavepacket['xFRF'], wavepacket['yFRF'], bathy, dx=1, dy=1, yBounds=[944, 947]) # non-inclusive index if you want 3 make 4 wide ## begin output @@ -108,10 +108,11 @@ def SwashSimSetup(startTime, inputDict): Dm=wavepacket['waveDm'], fileNameBase=date_str, path_prefix=path_prefix, version_prefix=version_prefix, nprocess=nprocess, runTime=runtime) # one compute core per cell in y + # write SWS file first swio.write_sws(swsinfo) swio.write_spec1D(wavepacket['freqbins'], wavepacket['fspec']) - swio.write_bot(gridDict['h']) + swio.write_bot(gridDict['elevation']) # now write QA/ swio.flags = None pickleName = os.path.join(path_prefix, date_str,'.pickle') @@ -123,7 +124,9 @@ def SwashAnalyze(startTime, inputDict, swio): """This runs the post process script for Swash wave will create plots and netcdf files at request Args: - inputDict (dict): this is an input dictionary that was generated with the + inputDict (dict): this is an input dictionary tvelopment +# Your branch is ahead of 'origin/development' by 8 commits. +# (use "git push" to publish your local commits)hat was generated with the keys from the project input yaml file startTime (str): input start time with datestring in format YYYY-mm-ddThh:mm:ssZ diff --git a/FUNWAVES_play.py b/potentialFilesForDeletion/FUNWAVES_play.py similarity index 100% rename from FUNWAVES_play.py rename to potentialFilesForDeletion/FUNWAVES_play.py diff --git a/funwavePlay.py b/potentialFilesForDeletion/funwavePlay.py similarity index 100% rename from funwavePlay.py rename to potentialFilesForDeletion/funwavePlay.py diff --git a/phaseResolvePlay.py b/potentialFilesForDeletion/phaseResolvePlay.py similarity index 100% rename from phaseResolvePlay.py rename to potentialFilesForDeletion/phaseResolvePlay.py diff --git a/prepdata b/prepdata index dbdbaf4..99bf957 160000 --- a/prepdata +++ b/prepdata @@ -1 +1 @@ -Subproject commit dbdbaf4d31b6718354ed43a4b31b4ae8c01b3463 +Subproject commit 99bf957a6168e6605db5e3b3da82afbaaa07488b