# -*- coding: utf-8 -*-
"""
Tools for standardizing the conversion of ADCP data to netCDF4 format for archiving and analysis.
"""
#@author: Chris Old
# IES, School of Engineering, The University of Edinburgh
#
#@note ... Tools developed under the RealTide project (WP2).
#
#@brief
# Tools for converting ReDAPT mat files to standardized netCDF.
#
#@details
# Contains functions for converting from ReDAPT RDI specific *.mat
# format files to standardized NetCDF4 format.
#
#@history 23/03/2021 -- Chris Old:
# Created original file on Tue Mar 23 16:25:01 2021
#-----------------------------------------------------------------------------
# IMPORTS
#-----------------------------------------------------------------------------
# Standard Python Dependencies
import os
import numpy as np
from datetime import datetime
# Non-Standard Python Dependencies
from netCDF4 import Dataset
import h5py
import utm
# Local Module Dependencies
from ocmw.core.timeFuncs import dateNumFromDateStr,matlab2python
from ocmw.core.fileTools import loadmat
# Other Dependencies
# NOTE:
# NetCDF4 in python requires the following packages installed
# cftime
# python_hdf4
# python_netcdf4
#
# These may come by default with an Anaconda install, otherwise the
# packages can be installed using pip.
#
#-----------------------------------------------------------------------------
# GLOBAL VARIABLES
#-----------------------------------------------------------------------------
ncfTimeUnits = 'days since 1900-01-01 00:00:00 UTC'
ncfEpochStr = '1900-01-01 00:00:00'
qc_tests = ['T6','T8','T9','T10','T14','T15','T18','T20','C1','C2']
qc_units = ['counts','counts','%','m/s','m/s','deg/s','counts','1/s','m','N/A']
qc_names = ['qc_signalStrength',
'qc_correlation',
'qc_percentGood',
'qc_velocityLimits',
'qc_errorVelocity',
'qc_rocTilts',
'qc_echoIntensity',
'qc_velocityGradient',
'qc_deployDepth',
'qc_internalRDIFlag']
qc_short_names = ['QC on signal strength',
'QC on correlation threshold',
'QC on percentage good threshold',
'QC on velocity limit thresholds',
'QC on error velocity',
'QC on rate of change of tilts',
'QC on echo intensity threshold',
'QC on velocity gradient threshold',
'QC on deployment depth threshold',
'QC on RDI internal bad data value']
qc_c2_names = ['qc_rdiFlag_east','qc_rdiFlag_north','qc_rdiFlag_vert',
'qc_rdiFlag_velx','qc_rdiFlag_vely','qc_rdiFlag_velz','qc_rdiFlag_vele',
'qc_rdiFlag_beam1','qc_rdiFlag_beam2','qc_rdiFlag_beam3','qc_rdiFlag_beam4']
qc_flag_values = [1, 2, 3, 4]
qc_flag_meaning = 'passed_qc_test not_evaluated questionable_quality failed_qc_test'
data_license = 'Creative Commons License: Attribution 4.0 International'
#-----------------------------------------------------------------------------
# CLASS DEFINITIONS
#-----------------------------------------------------------------------------
#-----------------------------------------------------------------------------
# FUNCTION DEFINITIONS
#-----------------------------------------------------------------------------
[docs]
def check_path(path):
pathExists = True
try:
os.makedirs(path)
except OSError:
if not os.path.isdir(path):
pathExists = False
raise
return pathExists
[docs]
def is_correct_file_extension(filename,extension):
fileName, fileExtension = os.path.splitext(filename)
correctFileExtension = False
if fileExtension == extension:
correctFileExtension = True
return correctFileExtension
[docs]
def set_file_extension(filename,extension):
fileName, fileExtension = os.path.splitext(filename)
newFileName = fileName + extension
return newFileName
[docs]
def check_file(path,filename,extension):
path = os.path.normpath(path)
check_path(path)
if not is_correct_file_extension(filename,extension):
filename = set_file_extension(filename,extension)
return path , filename
[docs]
def create_adcp_netcdf(ncfilename, num_beams=None, num_bins=None, num_recs=None):
"""
Create a generic netCDF data file structure
"""
# create new netcdf file
ncfile = Dataset(ncfilename,'w',format='NETCDF4')
# define dimensions
ncfile.createDimension('ndof', 3) # number of degrees of freedom for instrument orientation
ncfile.createDimension('nbeams', num_beams)
ncfile.createDimension('nbins', num_bins)
ncfile.createDimension('nrecs', num_recs)
# define variables
#--------------------------------
# scalars
#--------------------------------
ncfile.createVariable('configuration',np.uint8)
lat = ncfile.createVariable('latitude',np.float32)
lat.standard_name = 'latitude'
lat.long_name = 'latitude'
lat.crs = 'WGS84'
lat.units = 'degrees N'
lon = ncfile.createVariable('longitude',np.float32)
lon.standard_name = 'longitude'
lon.long_name = 'longitude'
lon.crs = 'WGS84'
lon.units = 'degrees E'
locx = ncfile.createVariable('easting',np.float32)
locx.standard_name = 'easting'
locx.long_name = 'easting'
locx.crs = 'UTM Zone 30N'
locx.units = 'm'
locy = ncfile.createVariable('northing',np.float32)
locy.standard_name = 'northing'
locy.long_name = 'northing'
locy.crs = 'UTM Zone 30N'
locy.units = 'm'
bmang = ncfile.createVariable('beam_angle',np.float32)
bmang.standard_name = 'beam angle'
bmang.long_name = 'beam angle from instrument vertical axis'
bmang.units = 'degrees'
num_beams = ncfile.createVariable('num_beams',np.int32)
num_beams.standard_name = 'number of beams'
num_beams.long_name = 'number of instrument transducers'
num_beams.units = ''
srate = ncfile.createVariable('sample_rate',np.float32)
srate.standard_name = 'sample rate'
srate.long_name = 'sampling rate'
srate.units = 'Hz'
headhab = ncfile.createVariable('head_hab',np.float32)
headhab.standard_name = 'head hab'
headhab.long_name = 'transducer head height above bed'
headhab.units = 'm'
blanking = ncfile.createVariable('blanking',np.float32)
blanking.standard_name = 'blanking distance'
blanking.long_name = 'blanking distance after transmit'
blanking.units = 'm'
bindist = ncfile.createVariable('bin_dist',np.float32)
bindist.standard_name = 'bin distance'
bindist.long_name = 'distance to bin 1'
bindist.units = 'm'
binsize = ncfile.createVariable('bin_size',np.float32)
binsize.standard_name = 'bin size'
binsize.long_name = 'profile bin size'
binsize.units = 'm'
num_bins = ncfile.createVariable('num_bins',np.int32)
num_bins.standard_name = 'number of bins'
num_bins.long_name = 'number of bins along profile'
num_bins.units = ''
#--------------------------------
# vectors
#--------------------------------
# timestamp
times = ncfile.createVariable('times',np.float64,('nrecs'))
times.standard_name = 'times'
times.long_name = 'record timestamp'
times.units = ncfTimeUnits
times.epoch_str = ncfEpochStr
times.calendar = 'gregorian'
# heading
heading = ncfile.createVariable('heading',np.float32,('nrecs'),fill_value=-99999.0)
heading.standard_name = 'heading'
heading.long_name = 'compass heading'
heading.units = 'degrees from north'
# pitch
pitch = ncfile.createVariable('pitch',np.float32,('nrecs'),fill_value=-99999.0)
pitch.standard_name = 'pitch'
pitch.long_name = 'tilt angle forward or aft'
pitch.units = 'degrees'
# roll
roll = ncfile.createVariable('roll',np.float32,('nrecs'),fill_value=-99999.0)
roll.standard_name = 'roll'
roll.long_name = 'tilt angle port or starboard'
roll.units = 'degrees'
# pressure
press = ncfile.createVariable('pressure',np.float32,('nrecs'),fill_value=-99999.0)
press.standard_name = 'pressure'
press.long_name = 'pressure at sensor'
press.units = 'dBar'
# temperature
temperature = ncfile.createVariable('temperature',np.float32,('nrecs'),fill_value=-99999.0)
temperature.standard_name = 'temperature'
temperature.long_name = 'tempearture at sensor'
temperature.units = 'degrees C'
# sound speed
sndspd = ncfile.createVariable('speedOfSound',np.float32,('nrecs'),fill_value=-99999.0)
sndspd.standard_name = 'sound speed'
sndspd.long_name = 'speed of sound in seawater'
sndspd.units = 'm/s'
# close output netcdf file
ncfile.close()
return
[docs]
def add_beam_amplitude(ncfile,num_beams):
"""
Add beam amplitude data structure
"""
for beam in range(num_beams):
var_name = 'ampl_beam'+str(beam+1)
ampl_beam = ncfile.createVariable(var_name,np.uint8,('nrecs','nbins'),fill_value=0)
ampl_beam.standard_name = 'beam amplitude'
ampl_beam.long_name = 'along-beam backscatter echo intensity'
ampl_beam.units = 'counts'
return
[docs]
def add_beam_correlation(ncfile,num_beams):
"""
Add beam correlation data structure
"""
for beam in range(num_beams):
var_name = 'corr_beam'+str(beam+1)
ampl_beam = ncfile.createVariable(var_name,np.uint8,('nrecs','nbins'),fill_value=0)
ampl_beam.standard_name = 'beam correlation'
ampl_beam.long_name = 'along-beam normalised echo auto-correlation at lag'
ampl_beam.units = 'counts'
return
[docs]
def add_beam_percent_good(ncfile,num_beams):
"""
Add beam percent good data structure
"""
for beam in range(num_beams):
var_name = 'pcgd_beam'+str(beam+1)
ampl_beam = ncfile.createVariable(var_name,np.uint8,('nrecs','nbins'),fill_value=0)
ampl_beam.standard_name = 'beam percent good'
ampl_beam.long_name = 'along-beam pergentage of good sub-pings'
ampl_beam.units = 'percentage'
return
[docs]
def add_beam_vel(ncfile,num_beams):
"""
Add beam velocity data structure
"""
for beam in range(num_beams):
var_name = 'vel_beam'+str(beam+1)
vel_beam = ncfile.createVariable(var_name,np.float32,('nrecs','nbins'),fill_value=-99999.0)
vel_beam.standard_name = 'beam velocity'
vel_beam.long_name = 'along-beam velocity'
vel_beam.sign_convention = 1.
vel_beam.units = 'm/s'
return
[docs]
def add_inst_vel(ncfile,num_beams=4):
"""
Add instrument coordinate velocity data structure
"""
vel_x = ncfile.createVariable('vel_instr_x',np.float32,('nrecs','nbins'),fill_value=-99999.0)
vel_x.standard_name = 'velocity x'
vel_x.long_name = 'velocity in instrument x coordinate direction'
vel_x.units = 'm/s'
vel_y = ncfile.createVariable('vel_instr_y',np.float32,('nrecs','nbins'),fill_value=-99999.0)
vel_y.standard_name = 'velocity y'
vel_y.long_name = 'velocity in instrument y coordinate direction'
vel_y.units = 'm/s'
if num_beams > 4:
vel_z1 = ncfile.createVariable('vel_instr_z1',np.float32,('nrecs','nbins'),fill_value=-99999.0)
vel_z1.standard_name = 'velocity z'
vel_z1.long_name = 'velocity in instrument z coordinate direction from transform'
vel_z1.units = 'm/s'
vel_z2 = ncfile.createVariable('vel_instr_z2',np.float32,('nrecs','nbins'),fill_value=-99999.0)
vel_z2.standard_name = 'velocity z'
vel_z2.long_name = 'velocity in instrument z coordinate direction from vertical beam'
vel_z2.units = 'm/s'
else:
vel_z = ncfile.createVariable('vel_instr_z',np.float32,('nrecs','nbins'),fill_value=-99999.0)
vel_z.standard_name = 'velocity z'
vel_z.long_name = 'velocity in instrument z coordinate direction'
vel_z.units = 'm/s'
vel_e = ncfile.createVariable('vel_instr_e',np.float32,('nrecs','nbins'),fill_value=-99999.0)
vel_e.standard_name = 'velocity error'
vel_e.long_name = 'error in instrument coordinate velocities'
vel_e.units = 'm/s'
return
[docs]
def add_earth_vel(ncfile):
"""
Add earth coordinate horizontal velocity data structure
"""
vel_east = ncfile.createVariable('vel_geo_east',np.float32,('nrecs','nbins'),fill_value=-99999.0)
vel_east.standard_name = 'eastward velocity'
vel_east.long_name = 'eastward component of velocity in geospatial Earth coordinates'
vel_east.units = 'm/s'
vel_north = ncfile.createVariable('vel_geo_north',np.float32,('nrecs','nbins'),fill_value=-99999.0)
vel_north.standard_name = 'northward velocity'
vel_north.long_name = 'northward component of velocity in geospatial Earth coordinates'
vel_north.units = 'm/s'
return
[docs]
def add_vertical_vel(ncfile):
"""
Add earth coordinate vertical velocity data structure
"""
vel_vert = ncfile.createVariable('vel_geo_vert',np.float32,('nrecs','nbins'),fill_value=-99999.0)
vel_vert.standard_name = 'vertical velocity'
vel_vert.long_name = 'vertical component of velocity in geospatial Earth coordinates'
vel_vert.units = 'm/s'
return
[docs]
def add_qc_data(ncfile):
"""
Add quality control data structure
"""
for indx, test in enumerate(qc_tests):
name = qc_names[indx]
dim2 = 'nbins'
if test == 'T15':
dim2 = 'ndof'
if test == 'C1':
ncfile.createVariable(name,np.uint8,('nrecs'),fill_value=0)
elif test == 'C2':
for name in qc_c2_names:
ncfile.createVariable(name,np.uint8,('nrecs',dim2),fill_value=0)
else:
ncfile.createVariable(name,np.uint8,('nrecs',dim2),fill_value=0)
[docs]
def add_stream_vel(ncfile):
"""
Add streamwise coordinate horizontal velocity data structure
"""
vel_stream = ncfile.createVariable('vel_stream',np.float32,('nrecs','nbins'),fill_value=-99999.0)
vel_stream.standard_name = 'streamwise velocity'
vel_stream.long_name = 'velocity parallel to principle flow direction'
vel_stream.units = 'm/s'
vel_trans = ncfile.createVariable('vel_trans',np.float32,('nrecs','nbins'),fill_value=-99999.0)
vel_trans.standard_name = 'transverse velocity'
vel_trans.long_name = 'velocity transverse to principle flow direction'
vel_trans.units = 'm/s'
return
[docs]
def add_adcp_array_fields(ncfilename, fields: list, num_beams):
"""
Add specific ADCP array field data structures
"""
# create new netcdf file
ncfile = Dataset(ncfilename,'r+',format='NETCDF4')
for field in fields:
#--------------------------------
# arrays
#--------------------------------
# beam velocities
if field == 'beam':
add_beam_vel(ncfile, num_beams)
# beam amplitude of backscatter
if field == 'ampl':
add_beam_amplitude(ncfile, num_beams)
# beam echo autocorrelation
if field == 'corr':
add_beam_correlation(ncfile, num_beams)
# beam percent good sub-pings
if field == 'pcgd':
add_beam_percent_good(ncfile, num_beams)
# instrument velocities
if field == 'inst':
add_inst_vel(ncfile,num_beams)
# earth velocities
if field == 'earth':
add_earth_vel(ncfile)
# vertical vel
if field == 'vert':
add_vertical_vel(ncfile)
# quality control data
if field == 'qc':
add_qc_data(ncfile)
# stream velocities
if field == 'strm':
add_stream_vel(ncfile)
# close output netcdf file
ncfile.close()
return
[docs]
def build_adcp_netcdf(path:str, filename: str, fields: list,
num_beams=None, num_bins=None, num_recs=None):
"""
Build a netCDF4 data file for a specific instrument
"""
# check path and file name
path , filename = check_file(path, filename, '.nc')
# create full file name for NetCDF file
ncfilename = os.path.join(path,filename)
create_adcp_netcdf(ncfilename, num_beams, num_bins, num_recs)
add_adcp_array_fields(ncfilename, fields, num_beams)
return
# Signature 5-beam instruments - pre-processing signature deployment software to *.mat
[docs]
def populate_sig_core(datpath,datfilename,ncpath,ncfilename,metadata):
"""
Populate the core data for a Nortek Signature NetCDF4 file from pre-processed data in matlab format.
"""
# create full file name for data file
data_file = os.path.join(datpath,datfilename).replace('\\','/')
# open source deployment data file
data = loadmat(data_file)
cfg = data['Config']
# check path and file name
ncpath , ncfilename = check_file(ncpath, ncfilename, '.nc')
# create full file name for NetCDF file
adcp_file = os.path.join(ncpath,ncfilename).replace('\\','/')
# open destination netcdf file
dst = Dataset(adcp_file,'r+',format='NETCDF4')
# set global attributes
# Description
dst.title = metadata.title
dst.projects = metadata.projects
dst.references = metadata.references
dst.site = metadata.site
dst.deploy_id= metadata.deploy_name
dst.geospatial_latitude = metadata.geospatial_latitude
dst.geospatial_longitude = metadata.geospatial_longitude
dst.geospatial_easting = metadata.geospatial_easting
dst.geospatial_northing = metadata.geospatial_northing
dst.geospatial_crs = metadata.geospatial_crs
dst.data_type = metadata.data_type
dst.processing_level = metadata.processing_level
dst.quality_control_indicator = metadata.quality_control_indicator
# Instrument
dst.instrument_type = metadata.instrument_type
dst.instrument_freq = metadata.instrument_freq
dst.instrument_serial_no = metadata.instrument_serial_no
# File generation
dst.contact = metadata.contact
dst.institution = metadata.institution
dst.data_assembly_center = metadata.data_assembly_center
if metadata.date_created == '':
dst.date_created = datetime.utcnow().strftime("%Y-%m-%d %H:%M:%S")
dst.date_updated = ''
else:
dst.date_updated = datetime.utcnow().strftime("%Y-%m-%d %H:%M:%S")
dst.history = metadata.history
dst.naming_authority = metadata.naming_authority
dst.data_format_version = metadata.data_format_version
dst.conventions = metadata.conventions
dst.netcdf_version = metadata.netcdf_version
if metadata.license == '':
dst.license = data_license
else:
dst.license = metadata.license
dst.citation = metadata.citation
dst.doi = metadata.doi
# Populate Configuration
dst.variables['configuration'].SerialNo = cfg['SerialNo']
dst.variables['configuration'].InstrumentName = cfg['InstrumentName']
dst.variables['configuration'].AnalogBoardRev = cfg['AnalogBoardRev']
dst.variables['configuration'].DigitalBoardRev = cfg['DigitalBoardRev']
dst.variables['configuration'].SensorBoardRev = cfg['SensorBoardRev']
dst.variables['configuration'].InterfaceBoardRev = cfg['InterfaceBoardRev']
dst.variables['configuration'].FPGAVersion = cfg['FPGAVersion']
dst.variables['configuration'].FWVersion = cfg['FWVersion']
dst.variables['configuration'].ComType = cfg['ComType']
dst.variables['configuration'].Baudrate = cfg['Baudrate']
dst.variables['configuration'].LedEnabled = cfg['LedEnabled']
dst.variables['configuration'].Orientation = cfg['Orientation']
dst.variables['configuration'].HeadFrequency = cfg['HeadFrequency']
dst.variables['configuration'].HasAccelerometer = np.uint8(eval(cfg['HasAccelerometer']))
dst.variables['configuration'].HasMagnetometer = np.uint8(eval(cfg['HasMagnetometer']))
dst.variables['configuration'].HasPressure = np.uint8(eval(cfg['HasPressure']))
dst.variables['configuration'].HasTemperature = np.uint8(eval(cfg['HasTemperature']))
dst.variables['configuration'].HasCompass = np.uint8(eval(cfg['HasCompass']))
dst.variables['configuration'].Hx = cfg['Hx']
dst.variables['configuration'].Hy = cfg['Hy']
dst.variables['configuration'].Hz = cfg['Hz']
dst.variables['configuration'].PressureOffset = cfg['PressureOffset']
dst.variables['configuration'].Declination = cfg['Declination']
dst.variables['configuration'].Plan_BurstEnabled = cfg['Plan_BurstEnabled']
dst.variables['configuration'].Plan_AverageEnabled = cfg['Plan_AverageEnabled']
dst.variables['configuration'].Plan_Salinity = cfg['Plan_Salinity']
dst.variables['configuration'].Plan_Filename = cfg['Plan_Filename']
dst.variables['configuration'].Plan_MaxVerticalVelocity = cfg['Plan_MaxVerticalVelocity']
dst.variables['configuration'].Plan_ProfileInterval = cfg['Plan_ProfileInterval']
dst.variables['configuration'].Plan_AverageDepthInterval = cfg['Plan_AverageDepthInterval']
dst.variables['configuration'].Plan_BurstInterval = cfg['Plan_BurstInterval']
dst.variables['configuration'].Plan_BurstDepthInterval = cfg['Plan_BurstDepthInterval']
dst.variables['configuration'].Plan_SoundVelocity = cfg['Plan_SoundVelocity']
dst.variables['configuration'].Plan_SerialOutputEnabled = cfg['Plan_SerialOutputEnabled']
dst.variables['configuration'].Burst_Altimeter = cfg['Burst_Altimeter']
dst.variables['configuration'].Burst_BlankingDistance = cfg['Burst_BlankingDistance']
dst.variables['configuration'].Burst_CellSize = cfg['Burst_CellSize']
dst.variables['configuration'].Burst_NBeams = cfg['Burst_NBeams']
dst.variables['configuration'].Burst_NCells = cfg['Burst_NCells']
dst.variables['configuration'].Burst_BeamSelection = cfg['Burst_BeamSelection']
dst.variables['configuration'].Burst_NSample = cfg['Burst_NSample']
dst.variables['configuration'].Burst_CoordSystem = cfg['Burst_CoordSystem']
dst.variables['configuration'].Burst_DataFormat = cfg['Burst_DataFormat']
dst.variables['configuration'].Burst_NPings = cfg['Burst_NPings']
dst.variables['configuration'].Burst_PowerLevel = cfg['Burst_PowerLevel']
dst.variables['configuration'].Burst_VelocityPrecision = cfg['Burst_VelocityPrecision']
dst.variables['configuration'].Burst_VelocityRange = cfg['Burst_VelocityRange']
dst.variables['configuration'].Burst_SamplingRate = cfg['Burst_SamplingRate']
dst.variables['configuration'].Burst_VelocityRangeBeam5 = cfg['Burst_VelocityRangeBeam5']
dst.variables['configuration'].Burst_RawAltimeter = cfg['Burst_RawAltimeter']
dst.variables['configuration'].Burst_HighResolution = np.uint8(eval(cfg['Burst_HighResolution']))
dst.variables['configuration'].Burst_HighResolution5 = np.uint8(eval(cfg['Burst_HighResolution5']))
dst.variables['configuration'].Burst_EchoSounder = np.uint8(eval(cfg['Burst_EchoSounder']))
dst.variables['configuration'].Burst_BottomTrack = np.uint8(eval(cfg['Burst_BottomTrack']))
dst.close()
return
[docs]
def populate_sig_fields(datpath,datfilename,ncpath,ncfilename,fields,indices,metadata):
"""
Populate the data fields for a Nortek Signature NetCDF4 file from pre-processed data in matlab format.
"""
# create full file name for data file
data_file = os.path.join(datpath,datfilename).replace('\\','/')
# open source deployment data file
data = loadmat(data_file)
cfg = data['Config']
src = data['Data']
# check path and file name
ncpath , ncfilename = check_file(ncpath, ncfilename, '.nc')
# create full file name for NetCDF file
adcp_file = os.path.join(ncpath,ncfilename).replace('\\','/')
# open destination netcdf file
dst = Dataset(adcp_file,'r+',format='NETCDF4')
startRec = dst.dimensions['nrecs'].size
# set global attributes
if hasattr(dst,'original_files'):
list(dst.original_files).append(datfilename)
list(dst.original_file_indices).append([indices[0],indices[-1]])
else:
dst.original_files = [datfilename]
dst.original_file_indices = [indices[0], indices[-1]]
# set variables
dst.variables['latitude'][0] = metadata.geospatial_latitude
dst.variables['longitude'][0] = metadata.geospatial_longitude
dst.variables['easting'][0] = metadata.geospatial_easting
dst.variables['northing'][0] = metadata.geospatial_northing
dst.variables['beam_angle'][0] = metadata.instrument_beam_angle
dst.variables['num_beams'][0] = cfg['Burst_NBeams']
dst.variables['sample_rate'][0] = cfg['Burst_SamplingRate']
dst.variables['blanking'][0] = cfg['Burst_BlankingDistance']
dst.variables['bin_dist'][0] = cfg['Burst_BlankingDistance']+cfg['Burst_CellSize']
dst.variables['bin_size'][0] = cfg['Burst_CellSize']
dst.variables['num_bins'][0] = cfg['Burst_NCells']
dst.variables['head_hab'][0] = metadata.instrument_head_hab
times = matlab2python(src['Burst_Time'][indices])
times = times - dateNumFromDateStr('1900-01-01 00:00:00')
dst.variables['times'][startRec:] = times
dst.variables['heading'][startRec:] = src['Burst_Heading'][indices]
dst.variables['pitch'][startRec:] = src['Burst_Pitch'][indices]
dst.variables['roll'][startRec:] = src['Burst_Roll'][indices]
dst.variables['pressure'][startRec:] = src['Burst_Pressure'][indices]
dst.variables['temperature'][startRec:] = src['Burst_Temperature'][indices]
dst.variables['speedOfSound'][startRec:] = src['Burst_Soundspeed'][indices]
for field in fields:
num_beams = np.int64(dst.variables['num_beams'][0])
if field == 'ampl':
for beam in range(num_beams):
if beam < 4:
varname = 'ampl_beam'+str(beam+1)
srcname = 'Burst_AmpBeam'+str(beam+1)
else:
varname = 'ampl_beam5'
srcname = 'IBurst_AmpBeam5'
dst.variables[varname][startRec:,:] = np.uint8(src[srcname][indices,:]*2.0)
if field == 'corr':
for beam in range(num_beams):
if beam < 4:
varname = 'corr_beam'+str(beam+1)
srcname = 'Burst_CorBeam'+str(beam+1)
else:
varname = 'corr_beam'+str(beam+1)
srcname = 'IBurst_CorBeam'+str(beam+1)
dst.variables[varname][startRec:,:] = np.uint8(src[srcname][indices,:])
if field == 'beam':
for beam in range(num_beams):
if beam < 4:
varname = 'vel_beam'+str(beam+1)
srcname = 'Burst_VelBeam'+str(beam+1)
else:
varname = 'vel_beam'+str(beam+1)
srcname = 'IBurst_VelBeam5'
dst.variables[varname][startRec:,:] = src[srcname][indices,:]
if field == 'inst':
dst.variables['vel_instr_x'][startRec:,:] = src['Burst_VelX'][indices,:]
dst.variables['vel_instr_y'][startRec:,:] = src['Burst_VelY'][indices,:]
vz1 = src['Burst_VelZ1'][indices,:]
vz2 = src['Burst_VelZ2'][indices,:]
if num_beams < 5:
dst.variables['vel_instr_z'][startRec:,:] = (vz1+vz2)*0.5
else:
dst.variables['vel_instr_z1'][startRec:,:] = (vz1+vz2)*0.5
dst.variables['vel_instr_z2'][startRec:,:] = src['IBurst_VelBeam5'][indices,:]
dst.variables['vel_instr_e'][startRec:,:] = np.abs(vz1-vz2)/np.sqrt(2.0)
if field == 'earth':
dst.variables['vel_geo_east'][startRec:,:] = src['Burst_VelEast'][indices,:]
dst.variables['vel_geo_north'][startRec:,:] = src['Burst_VelNorth'][indices,:]
vu1 = src['Burst_VelUp1'][indices,:]
vu2 = src['Burst_VelUp2'][indices,:]
dst.variables['vel_geo_vert'][startRec:,:] = (vu1+vu2)*0.5
dst.close()
return
# RDI Workhorse - ReDAPT pre-processing to *.mat files
[docs]
def populate_adcp_core(datpath,datfilename,ncpath,ncfilename):
"""
Populate the core data for an RDI Workhorse NetCDF4 file from pre-processed data in ReDAPT matlab format.
"""
# create full file name for data file
data_file = os.path.join(datpath,datfilename).replace('\\','/')
# open source deployment data file
data = h5py.File(data_file, 'r')
src = data['data']
# check path and file name
ncpath , ncfilename = check_file(ncpath, ncfilename, '.nc')
# create full file name for NetCDF file
adcp_file = os.path.join(ncpath,ncfilename).replace('\\','/')
# open destination netcdf file
dst = Dataset(adcp_file,'r+',format='NETCDF4')
deploy_name = datpath.replace('\\','/').split('/')[-1]
# set global attributes
# Description
dst.title = 'ReDAPT RDI ADCP data standardization'
dst.projects = 'ReDAPT - Reliable Data Acquisition Platform for Tidal'
dst.references = 'https://www.realtide.eu/ http://redapt.eng.ed.ac.uk/'
dst.site = 'EMEC Test Site, Falls of Warness, Orkney'
dst.deploy_id= deploy_name
dst.geospatial_easting = src['Config']['Location_East'][0][0]
dst.geospatial_northing = src['Config']['Location_North'][0][0]
dst.geospatial_crs = 'WGS 84 / '+''.join(chr(i) for i in (np.transpose(src['Config']['Location_Reference'][:])[0]))
dst.data_type = 'field measurement'
dst.processing_level = '1'
dst.quality_control_indicator = '1'
# Instrument
dst.instrument_type = 'RDI Workhorse'
dst.instrument_freq = '600kHz'
dst.instrument_serial_no = np.int32(src['Config']['instSerialNumber'][0][0]).astype(str)
# File generation
dst.contact = ''
dst.institution = 'IES, School of Engineering, The University of Edinburgh, UK'
dst.data_assembly_center = ''
dst.date_created = datetime.utcnow().strftime("%Y-%m-%d %H:%M:%S")
dst.date_updated = ''
dst.history = ''
dst.naming_authority = 'University of Edinburgh'
dst.data_format_version = '1.0'
dst.conventions = 'CF-1.8'
dst.netcdf_version = 'netCDF-4 classic model'
dst.license = data_license
dst.citation = ''
dst.doi = ''
# set configuration
dst.variables['configuration'].beam_config = 'divergent'
dst.variables['configuration'].beam_angle = np.int32(src['Config']['beamAngle'][0][0])
dst.variables['configuration'].num_beams = np.int32(src['Config']['numBeams'][0][0])
tppm = np.float32(src['Config']['tppMinutes'][0][0])*60.0
tpps = np.float32(src['Config']['tppSeconds'][0][0])
tpph = np.float32(src['Config']['tppHundredths'][0][0])/100.0
ppens = np.float32(src['Config']['pingsPerEnsemble'][0][0])
sample_rate = 1.0/(ppens*(tppm+tpps+tpph))
dst.variables['configuration'].sample_rate = sample_rate
dst.variables['configuration'].bin_dist = np.int32(src['Config']['bin1Distance'][0][0])
dst.variables['configuration'].bin_size = np.int32(src['Config']['depthCellLength'][0][0])
dst.variables['configuration'].num_bins = np.int32(src['Config']['numCells'][0][0])
dst.variables['configuration'].blanking = np.int32(src['Config']['blankAfterTransmit'][0][0])
dst.variables['configuration'].sensorsAvailable = ''.join(str(i) for i in (np.transpose(src['Config']['sensorsAvailable'][:])[0]))
dst.variables['configuration'].sensorSource = ''.join(str(i) for i in (np.transpose(src['Config']['sensorSource'][:])[0]))
dst.variables['configuration'].systemConfiguration = ''.join(str(i) for i in (np.transpose(src['Config']['systemConfiguration'][:])[0]))
dst.variables['configuration'].systemPower = np.int32(src['Config']['systemPower'][0][0])
dst.variables['configuration'].systemBandwidth = np.int32(src['Config']['systemBandwidth'][0][0])
# Sampling
dst.variables['configuration'].profilingMode = np.int32(src['Config']['profilingMode'][0][0])
dst.variables['configuration'].pingPerEnsemble = np.int32(src['Config']['pingsPerEnsemble'][0][0])
dst.variables['configuration'].tppMinutes = np.int32(src['Config']['tppMinutes'][0][0])
dst.variables['configuration'].tppSeconds = np.int32(src['Config']['tppSeconds'][0][0])
dst.variables['configuration'].tppHundredths = np.int32(src['Config']['tppHundredths'][0][0])
# Internal processing
dst.variables['configuration'].fixedLeaderId = np.int32(src['Config']['fixedLeaderId'][0][0])
dst.variables['configuration'].coordinateTransform = ''.join(str(i) for i in (np.transpose(src['Config']['coordinateTransform'][:])[0]))
dst.variables['configuration'].headingAlignment = np.float32(src['Config']['headingAlignment'][0][0])
dst.variables['configuration'].headingBias = np.float32(src['Config']['headingBias'][0][0])
dst.variables['configuration'].errVelocityMax = np.int32(src['Config']['errVelocityMax'][0][0])
dst.variables['configuration'].falseTargetThresh = np.int32(src['Config']['falseTargetThresh'][0][0])
dst.variables['configuration'].lowCorrThresh = np.int32(src['Config']['lowCorrThresh'][0][0])
dst.variables['configuration'].gdMinimum = np.int32(src['Config']['gdMinimum'][0][0])
dst.variables['configuration'].wpRefLayerAvgStart = np.int32(src['Config']['wpRefLayerAvgStart'][0][0])
dst.variables['configuration'].wpRefLayerAvgEnd = np.int32(src['Config']['wpRefLayerAvgEnd'][0][0])
dst.variables['configuration'].realSimFlag = np.int32(src['Config']['realSimFlag'][0][0])
# internal settings
dst.variables['configuration'].numCodeReps = np.int32(src['Config']['numCodeReps'][0][0])
dst.variables['configuration'].lagLength = np.int32(src['Config']['lagLength'][0][0])
dst.variables['configuration'].transmitLagDistance = np.int32(src['Config']['transmitLagDistance'][0][0])
dst.variables['configuration'].xmitPulseLength = np.int32(src['Config']['xmitPulseLength'][0][0])
# set serial numbers etc.
dst.variables['configuration'].serialNumber = np.int32(src['Config']['instSerialNumber'][0][0]).astype(str)
dst.variables['configuration'].cpuBoardSerialNo = np.int64(src['Config']['cpuBoardSerialNo'][0][0]).astype(str)
dst.variables['configuration'].cpuFirmwareVersion = np.int32(src['Config']['cpuFirmwareVersion'][0][0]).astype(str)
dst.variables['configuration'].cpuFirmwareRevision = np.int32(src['Config']['cpuFirmwareRevision'][0][0]).astype(str)
dst.close()
data.close()
return
[docs]
def populate_adcp_fields(datpath,datfilename,ncpath,ncfilename,fields,indices):
"""
Populate the data fields for an RDI Workhorse NetCDF4 file from pre-processed data in ReDAPT matlab format.
"""
# create full file name for data file
data_file = os.path.join(datpath,datfilename).replace('\\','/')
# open source deployment data file
data = h5py.File(data_file, 'r')
src = data['data']
# check path and file name
ncpath , ncfilename = check_file(ncpath, ncfilename, '.nc')
# create full file name for NetCDF file
adcp_file = os.path.join(ncpath,ncfilename).replace('\\','/')
# open destination netcdf file
dst = Dataset(adcp_file,'r+',format='NETCDF4')
startRec = dst.dimensions['nrecs'].size
# set global attributes
if hasattr(dst,'original_files'):
list(dst.original_files).append(datfilename)
list(dst.original_file_indices).append([indices[0],indices[-1]])
else:
dst.original_files = [datfilename]
dst.original_file_indices = [indices[0], indices[-1]]
# set variables
dst.variables['easting'][0] = src['Config']['Location_East'][0][0]
dst.variables['northing'][0] = src['Config']['Location_North'][0][0]
dst.variables['beam_angle'][0] = src['Config']['beamAngle'][0][0]
dst.variables['num_beams'][0] = np.int32(src['Config']['numBeams'][0][0])
dst.variables['sample_rate'][0] = np.round(src['Config']['Sample_Rate'][0][0]*10000.0)/10000.0
dst.variables['blanking'].blanking = src['Config']['blankAfterTransmit'][0][0]/100.0
dst.variables['bin_dist'][0] = src['Config']['bin1Distance'][0][0]/100.0
dst.variables['bin_size'][0] = src['Config']['depthCellLength'][0][0]/100.0
dst.variables['num_bins'][0] = np.int32(src['Config']['numCells'][0][0])
dst.variables['head_hab'][0] = 0.8
times = matlab2python(np.squeeze(src['TimestampUTC'])[indices])
times = times - dateNumFromDateStr('1900-01-01 00:00:00')
dst.variables['times'][startRec:] = times
dst.variables['heading'][startRec:] = np.squeeze(src['Heading'])[indices]
dst.variables['pitch'][startRec:] = np.squeeze(src['Pitch'])[indices]
dst.variables['roll'][startRec:] = np.squeeze(src['Roll'])[indices]
dst.variables['pressure'][startRec:] = np.squeeze(src['Pressure'])[indices]
dst.variables['temperature'][startRec:] = np.squeeze(src['Temperature'])[indices]
dst.variables['speedOfSound'][startRec:] = np.squeeze(src['SpeedOfSound'])[indices]
for field in fields:
num_beams = np.int64(dst.variables['num_beams'][0])
if field == 'ampl':
for beam in range(num_beams):
varname = 'ampl_beam'+str(beam+1)
srcname = 'Amp_Beam'+str(beam+1)
dst.variables[varname][startRec:,:] = np.uint8(src[srcname][:,indices].transpose())
if field == 'corr':
for beam in range(num_beams):
varname = 'corr_beam'+str(beam+1)
srcname = 'Corr_Beam'+str(beam+1)
dst.variables[varname][startRec:,:] = np.uint8(src[srcname][:,indices].transpose())
if field == 'pcgd':
for beam in range(num_beams):
varname = 'pcgd_beam'+str(beam+1)
srcname = 'PercentGood_Beam'+str(beam+1)
dst.variables[varname][startRec:,:] = np.uint8(src[srcname][:,indices].transpose())
if field == 'beam':
for beam in range(num_beams):
varname = 'vel_beam'+str(beam+1)
srcname = 'Vel_Beam'+str(beam+1)
dst.variables[varname][startRec:,:] = src[srcname][:,indices].transpose()
if field == 'inst':
dst.variables['vel_instr_x'][startRec:,:] = src['Vel_Instr_X'][:,indices].transpose()
dst.variables['vel_instr_y'][startRec:,:] = src['Vel_Instr_Y'][:,indices].transpose()
dst.variables['vel_instr_z'][startRec:,:] = src['Vel_Instr_Z'][:,indices].transpose()
dst.variables['vel_instr_e'][startRec:,:] = src['Vel_Instr_Error'][:,indices].transpose()
if field == 'earth':
dst.variables['vel_geo_east'][startRec:,:] = src['Vel_Geo_East'][:,indices].transpose()
dst.variables['vel_geo_north'][startRec:,:] = src['Vel_Geo_North'][:,indices].transpose()
dst.variables['vel_geo_vert'][startRec:,:] = src['Vel_Geo_Up'][:,indices].transpose()
if field == 'strm':
dst.variables['vel_stream'][startRec:,:] = src['U'][:,indices].transpose()
dst.variables['vel_trans'][startRec:,:] = src['V'][:,indices].transpose()
if field == 'vert':
dst.variables['vel_vert'][startRec:,:] = src['W'][:,indices].transpose()
if field == 'qc':
for indx, test in enumerate(qc_tests):
name = qc_names[indx]
std_name = qc_short_names[indx]
long_name = ''.join(chr(k) for k in np.squeeze(src['QCconfig'][test]['name']))
info = ''.join(chr(k) for k in np.squeeze(src['QCconfig'][test]['info']))
units = qc_units[indx]
if test == 'C1':
dst.variables[name][startRec:] = src['QC'][test][:,indices].transpose()
if not hasattr(dst.variables[name],'standard_name'):
dst.variables[name].standard_name = std_name
dst.variables[name].long_name = long_name
dst.variables[name].description = info
dst.variables[name].thresh_min = np.float32(src['QCconfig'][test]['thresh_min'][0][0])
dst.variables[name].thresh_max = np.float32(src['QCconfig'][test]['thresh_max'][0][0])
dst.variables[name].units = units
dst.variables[name].flag_values = qc_flag_values
dst.variables[name].flag_meanings = qc_flag_meaning
elif test == 'C2':
dst.variables['qc_rdiFlag_east'][startRec:,:] = src['QC'][test]['Vel_Geo_East'][:,indices].transpose()
dst.variables['qc_rdiFlag_north'][startRec:,:] = src['QC'][test]['Vel_Geo_North'][:,indices].transpose()
dst.variables['qc_rdiFlag_vert'][startRec:,:] = src['QC'][test]['Vel_Geo_Up'][:,indices].transpose()
dst.variables['qc_rdiFlag_velx'][startRec:,:] = src['QC'][test]['Vel_Instr_X'][:,indices].transpose()
dst.variables['qc_rdiFlag_vely'][startRec:,:] = src['QC'][test]['Vel_Instr_Y'][:,indices].transpose()
dst.variables['qc_rdiFlag_velz'][startRec:,:] = src['QC'][test]['Vel_Instr_Z'][:,indices].transpose()
dst.variables['qc_rdiFlag_vele'][startRec:,:] = src['QC'][test]['Vel_Instr_Error'][:,indices].transpose()
dst.variables['qc_rdiFlag_beam1'][startRec:,:] = src['QC'][test]['Vel_Beam1'][:,indices].transpose()
dst.variables['qc_rdiFlag_beam2'][startRec:,:] = src['QC'][test]['Vel_Beam2'][:,indices].transpose()
dst.variables['qc_rdiFlag_beam3'][startRec:,:] = src['QC'][test]['Vel_Beam3'][:,indices].transpose()
dst.variables['qc_rdiFlag_beam4'][startRec:,:] = src['QC'][test]['Vel_Beam4'][:,indices].transpose()
if not hasattr(dst.variables['qc_rdiFlag_east'],'standard_name'):
dst.variables['qc_rdiFlag_east'].standard_name = std_name
dst.variables['qc_rdiFlag_east'].long_name = long_name
dst.variables['qc_rdiFlag_east'].description = info
dst.variables['qc_rdiFlag_east'].units = units
dst.variables['qc_rdiFlag_east'].flag_values = qc_flag_values
dst.variables['qc_rdiFlag_east'].flag_meanings = qc_flag_meaning
dst.variables['qc_rdiFlag_north'].standard_name = std_name
dst.variables['qc_rdiFlag_north'].long_name = long_name
dst.variables['qc_rdiFlag_north'].description = info
dst.variables['qc_rdiFlag_north'].units = units
dst.variables['qc_rdiFlag_north'].flag_values = qc_flag_values
dst.variables['qc_rdiFlag_north'].flag_meanings = qc_flag_meaning
dst.variables['qc_rdiFlag_vert'].standard_name = std_name
dst.variables['qc_rdiFlag_vert'].long_name = long_name
dst.variables['qc_rdiFlag_vert'].description = info
dst.variables['qc_rdiFlag_vert'].units = units
dst.variables['qc_rdiFlag_vert'].flag_values = qc_flag_values
dst.variables['qc_rdiFlag_vert'].flag_meanings = qc_flag_meaning
dst.variables['qc_rdiFlag_velx'].standard_name = std_name
dst.variables['qc_rdiFlag_velx'].long_name = long_name
dst.variables['qc_rdiFlag_velx'].description = info
dst.variables['qc_rdiFlag_velx'].units = units
dst.variables['qc_rdiFlag_velx'].flag_values = qc_flag_values
dst.variables['qc_rdiFlag_velx'].flag_meanings = qc_flag_meaning
dst.variables['qc_rdiFlag_vely'].standard_name = std_name
dst.variables['qc_rdiFlag_vely'].long_name = long_name
dst.variables['qc_rdiFlag_vely'].description = info
dst.variables['qc_rdiFlag_vely'].units = units
dst.variables['qc_rdiFlag_vely'].flag_values = qc_flag_values
dst.variables['qc_rdiFlag_vely'].flag_meanings = qc_flag_meaning
dst.variables['qc_rdiFlag_velz'].standard_name = std_name
dst.variables['qc_rdiFlag_velz'].long_name = long_name
dst.variables['qc_rdiFlag_velz'].description = info
dst.variables['qc_rdiFlag_velz'].units = units
dst.variables['qc_rdiFlag_velz'].flag_values = qc_flag_values
dst.variables['qc_rdiFlag_velz'].flag_meanings = qc_flag_meaning
dst.variables['qc_rdiFlag_vele'].standard_name = std_name
dst.variables['qc_rdiFlag_vele'].long_name = long_name
dst.variables['qc_rdiFlag_vele'].description = info
dst.variables['qc_rdiFlag_vele'].units = units
dst.variables['qc_rdiFlag_vele'].flag_values = qc_flag_values
dst.variables['qc_rdiFlag_vele'].flag_meanings = qc_flag_meaning
dst.variables['qc_rdiFlag_beam1'].standard_name = std_name
dst.variables['qc_rdiFlag_beam1'].long_name = long_name
dst.variables['qc_rdiFlag_beam1'].description = info
dst.variables['qc_rdiFlag_beam1'].units = units
dst.variables['qc_rdiFlag_beam1'].flag_values = qc_flag_values
dst.variables['qc_rdiFlag_beam1'].flag_meanings = qc_flag_meaning
dst.variables['qc_rdiFlag_beam2'].standard_name = std_name
dst.variables['qc_rdiFlag_beam2'].long_name = long_name
dst.variables['qc_rdiFlag_beam2'].description = info
dst.variables['qc_rdiFlag_beam2'].units = units
dst.variables['qc_rdiFlag_beam2'].flag_values = qc_flag_values
dst.variables['qc_rdiFlag_beam2'].flag_meanings = qc_flag_meaning
dst.variables['qc_rdiFlag_beam3'].standard_name = std_name
dst.variables['qc_rdiFlag_beam3'].long_name = long_name
dst.variables['qc_rdiFlag_beam3'].description = info
dst.variables['qc_rdiFlag_beam3'].units = units
dst.variables['qc_rdiFlag_beam3'].flag_values = qc_flag_values
dst.variables['qc_rdiFlag_beam3'].flag_meanings = qc_flag_meaning
dst.variables['qc_rdiFlag_beam4'].standard_name = std_name
dst.variables['qc_rdiFlag_beam4'].long_name = long_name
dst.variables['qc_rdiFlag_beam4'].description = info
dst.variables['qc_rdiFlag_beam4'].units = units
dst.variables['qc_rdiFlag_beam4'].flag_values = qc_flag_values
dst.variables['qc_rdiFlag_beam4'].flag_meanings = qc_flag_meaning
else:
dst.variables[name][startRec:,:] = src['QC'][test][:,indices].transpose()
if not hasattr(dst.variables[name],'standard_name'):
dst.variables[name].standard_name = std_name
dst.variables[name].long_name = long_name
dst.variables[name].description = info
if test in ['T6','T8','T9']:
dst.variables[name].threshold = np.float32(src['QCconfig'][test]['thresh_min'][0][0])
elif test in ['T18','T20']:
dst.variables[name].threshold = np.float32(src['QCconfig'][test]['thresh_max'][0][0])
else:
dst.variables[name].thresh_min = np.float32(src['QCconfig'][test]['thresh_min'][0][0])
dst.variables[name].thresh_max = np.float32(src['QCconfig'][test]['thresh_max'][0][0])
dst.variables[name].units = units
dst.variables[name].flag_values = qc_flag_values
dst.variables[name].flag_meanings = qc_flag_meaning
dst.close()
return
#-----------------------------------------------------------------------------
# MODULE TEST CODE
#-----------------------------------------------------------------------------