#!/usr/bin/python #======================================================================= Scope = 'Read in data radial data files and average consecutive gates.\n'+\ 'compresses to byte format for udp packets and write out to\n'+\ 'file for sending.' # # This code is owned by Remote Sensing Solutions. It is intended for use # by Dr. Paul Chang of NOAA/NESDIS to provide NCO with real-time radial data # from the NOAA WP-3D Tail Doppler Radar. It is to be deployed on the NOAA # Aircraft Operations server. Use beyond the above purpose or distribution # beyond AOC requires approval from Dr. James Carswell of Remote Sensing # Solutions, 3179 Main Street, Barnstable, MA 02630. # # Written by: \ Author = 'James Carswell and Dan Robinson' Institution = 'Remote Sensing Solutions, Inc.' # Version versionTag = '0.1.1' # # Modifications # (List version, author, date and then changes.) # # V0.1.1 - Jim Carswell, 20100208 # 1) Added check when calculating numberOfAveBins to use the smaller of # numberOfBins (in header) and maxNumberOfGates. # #======================================================================= # Import Modules import time import string import os.path import dircache import struct import array import sys import zlib import gzip import zipfile def openNewFile(path, missionId, aircraftId, TDNumber, myTime): strAircraft = '?' if aircraftId == 43: strAircraft = 'I' elif aircraftId == 42: strAircraft = 'H' elif aircraftId == 49: strAircraft = 'N' fileName = '%08d%s%d-%s.dat' % \ (missionId,strAircraft,TDNumber, time.strftime('%H%M', myTime)) print '\nCreating file ' + fileName outputFile = open(os.path.join(path, fileName),'wb') return outputFile, fileName, strAircraft def addToFilesWritten(logPath, fileName): # Check for Output data file path if not(os.path.exists(logPath)): print '-> Log file path: '+logPath+' does not exist. Creating it.' os.makedirs(logPath) f = open(os.path.join(logPath, "tailRadarFilesWritten.list"), "a+") f.write(fileName + "\n") f.close() def processRadialFile( radialFilePath, \ newFileFlag, \ outputFilePath, \ archivePath, \ numberOfGatesToAverage, \ maxNumberOfGates, \ filesToCompress, \ verbose): headerTag = 0 dataTag = 1 stdOutCnt = " " radialFile = open(radialFilePath,mode='rb') if (radialFilePath.endswith('.gz')): # Open Underlying File zippedFile = open(radialFilePath,mode='rb') # Create the uncompressing wrapper for the zipped file # This allow us to treat the file as if it was not zipped radialFile = gzip.GzipFile(fileobj=zippedFile, mode='rt') # Read Input File Header header='' for lineNumber in xrange(8): header=header+radialFile.readline() splitHeader = header.split('\n') splitHeader = splitHeader[0:8] # Build Output Header Structure missionId = long(splitHeader[0]) aircraftId = int(splitHeader[1]) TDNumber = int(splitHeader[2]) numberOfBins = int(splitHeader[3]) if numberOfBins <= maxNumberOfGates: numberOfAveBins = int(int(splitHeader[3])/numberOfGatesToAverage) else: numberOfAveBins = int(int(maxNumberOfGates)/numberOfGatesToAverage) rangeDelay = float(splitHeader[4]) rangeGateResolution = float(splitHeader[5])*numberOfGatesToAverage reflectivityFlag = int(splitHeader[6][0]) DopplerFlag = int(splitHeader[6][1]) spectralWidthFlag = int(splitHeader[6][2]) dataType = reflectivityFlag + DopplerFlag*2 + \ spectralWidthFlag*4 antennaNumber = int(splitHeader[7]) headerFormat = '=l3h2f2h' packedHeader = struct.pack(headerFormat, missionId, aircraftId, TDNumber, numberOfAveBins, \ rangeDelay, rangeGateResolution, dataType, antennaNumber ) # Determine data block read size in bytes numberOfParameters = reflectivityFlag + DopplerFlag + spectralWidthFlag numberOfBinsPerLine = 10 numberOfCharactersPerBin = 6 numberOfLineFeeds = int(numberOfBins/numberOfBinsPerLine) if numberOfBins % numberOfBinsPerLine <> 0: numberOfLineFeeds = numberOfLineFeeds+1 numberOfBytesPerBlock = numberOfBins*numberOfCharactersPerBin + \ numberOfLineFeeds eof = 0 blockCount = long(-1) headerBlockCount = 0 rHeader = radialFile.readline() bytesWritten = 0 lastTime = time.time() startTime = lastTime z = 0 outputBuffer = array.array('B') tlvStructHeader = '=2HL' tlvStructHeaderSize = struct.calcsize(tlvStructHeader) minPayload = 1500 maxPayload = 0 minCount = 0 maxCount = 0 payloadLen = 0 readCount = 0 firstPointTime = 0.0 newFileTime = 300.0 # 5 minutes fileCreated = 0 bulletinPrefix = '****0000000000****\n' prefixLengthOffset = 4 # start of length field in prefix fileCount = 0 compFileStart = 0 fileList = [] zipFileNameFormat = '%d%s%d_TA_%d_%s-%s.zip' # Loop through the file while rHeader: radarHeaderFormat = '=3hd4f' radarHeaderLen = struct.calcsize(radarHeaderFormat) # Read in data block radarHeader = rHeader.split() unixTime = float(radarHeader[0]) latitude = float(radarHeader[1]) longitude = float(radarHeader[2]) altitude = int(float(radarHeader[3])) azimuth = float(radarHeader[4]) incidence = float(radarHeader[5]) # Read in radar data block readBlock = radialFile.read(numberOfBytesPerBlock) if (readBlock): dataBlock = readBlock.split() # Process Doppler Data if DopplerFlag: Doppler = array.array('f') if len(dataBlock) < numberOfBins+1: # Convert to float array n = len(dataBlock) if n > maxNumberOfGates: n = maxNumberOfGates for bin in xrange(n): Doppler = Doppler + \ array.array('f',[float(dataBlock[bin].strip('\n'))]) # Average data and create output array (short, tenths of m/s) aveDoppler = array.array('h') for aveBin in xrange(numberOfAveBins): tmp=0. n=0. for i in xrange(numberOfGatesToAverage): bin = (aveBin*numberOfGatesToAverage)+i if abs(Doppler[bin]) < 4096: n=n+1 tmp = tmp+Doppler[bin] if n > 0: tmp1 = int((tmp/n+4096.)+.5) aveDoppler = aveDoppler + \ array.array('h', [int(tmp1-4096.)]) else: aveDoppler = aveDoppler + \ array.array('h',[-8888]) # Write radar header radarHeaderOut = struct.pack('=dfffff', unixTime, latitude, longitude, altitude,\ azimuth, incidence) # Create new file? timeDiff = unixTime - firstPointTime if timeDiff > newFileTime or fileCreated == 0: if fileCreated == 1: # if at end of file, write length to header, # close and send to destination (separate script) outputFile.flush() # include length of URNT16 header, but not prefix totalLength = bytesWritten + len(urnt16Header) outputFile.seek( prefixLengthOffset ) strLength = '%010d' % (totalLength) outputFile.write( strLength ) # move pointer to EOF and close outputFile.seek( 0,2 ) print '\nFile complete. Closing file.' outputFile.close( ) addToFilesWritten( os.path.join(outputFilePath, 'log'), fname ) fileList.append( fname ) # Handle file compression fileCount += 1 print 'file count', fileCount if filesToCompress == 1 or fileCount % filesToCompress == 1: print 'first file' compFileStart = time.strftime('%H%M', myTime) if fileCount % filesToCompress == 0: # Do compression zipFileName = zipFileNameFormat % (missionId, strAircraft, TDNumber, filesToCompress, compFileStart, time.strftime('%H%M', time.gmtime(unixTime))) print 'Creating archive:', zipFileName zipFilePath = os.path.join(archivePath, zipFileName) arcfile = zipfile.ZipFile(zipFilePath, "w") for i in range(len(fileList)): arcfile.write(os.path.join(outputFilePath, fileList[i]), fileList[i], zipfile.ZIP_DEFLATED) arcfile.close() fileList = [] myTime = time.gmtime(unixTime) outputFile, fname, strAircraft= openNewFile( outputFilePath, missionId, aircraftId, TDNumber, myTime ) outputFile.write( bulletinPrefix ) # Append DateTimeGroup to URNT16 Header # Use time of first data point tupleTime = time.gmtime(unixTime) urnt16Header = 'URNT16 KWBC ' + time.strftime('%d%H%M', tupleTime) + '\n' outputFile.write( urnt16Header ) outputFile.write( packedHeader ) bytesWritten = len( packedHeader ) startTime = time.time() firstPointTime = unixTime firstDataPoint = 0 fileCreated = 1 bytesWritten += len(radarHeaderOut)+(len(aveDoppler)*aveDoppler.itemsize) outputFile.write( radarHeaderOut ) outputFile.write( aveDoppler ) outputFile.flush() # Block Counter if verbose ==1: if blockCount % 10e3 == 0: thisTime = time.time() print 'Block Count:',blockCount, ' Processing Time: ',thisTime-lastTime lastTime = thisTime else: if blockCount % 1e3 == 0: nbs = len(stdOutCnt) stdOutCnt = 'Block Count: %s' % blockCount sys.stdout.write('\b'*nbs+stdOutCnt) sys.stdout.flush() else: print '-> Bad Radar Block %s'%blockCount print len(dataBlock), numberOfBins+1 # Read the radar header, if any rHeader = radialFile.readline() # Done looping through file, clean up currently open file and last archive if fileCreated == 1: # if at end of file, write length to header, # close and send to destination (separate script) outputFile.flush() # include length of URNT16 header, but not prefix totalLength = bytesWritten + len(urnt16Header) outputFile.seek( prefixLengthOffset ) strLength = '%010d' % (totalLength) outputFile.write( strLength ) # move pointer to EOF and close outputFile.seek( 0,2 ) print '\nFile complete. Closing file.' outputFile.close( ) addToFilesWritten( os.path.join(outputFilePath, 'log'), fname ) fileList.append( fname ) # Handle file compression fileCount += 1 print 'file count', fileCount if filesToCompress == 1 or fileCount % filesToCompress == 1: print 'first file' compFileStart = time.strftime('%H%M', myTime) # Do compression if fileCount % filesToCompress == 0: compFileCount = filesToCompress else: compFileCount = fileCount % filesToCompress zipFileName = zipFileNameFormat % (missionId, strAircraft, TDNumber, compFileCount, compFileStart, time.strftime('%H%M', time.gmtime(unixTime))) print 'Creating archive:', zipFileName zipFilePath = os.path.join(archivePath, zipFileName) arcfile = zipfile.ZipFile(zipFilePath, "w") for i in range(len(fileList)): arcfile.write(os.path.join(outputFilePath, fileList[i]), fileList[i], zipfile.ZIP_DEFLATED) arcfile.close()