root/DPWP/trunk/DPWP/ADCP_splitter/adcpwaves.py

"""This code takes the raw waves portion of the adcp data from the splitter
and outputs pressure, range, orbital velocity and system info text files that
can be read into matlab.

NOTE:  Transducer height above bottom is hardwired into this code at .4m.
Need to update code in the future to be a user selected option"""




import sys, struct, math, shutil, os

try:
    infileName = sys.argv[1]

except:
    print 'error'
    sys.exit(1)


ifile = open(infileName, 'rb')

# set the transducer height above bottom
thab = 0.4


#header
firstHeader=31103
Header=firstHeader
while Header == firstHeader:
    readFirstHeader=ifile.read(2)
    if len(readFirstHeader) == 0:
        print('\n End of file')
        break

    pressureOfile = open('pressure', 'w')
    rangeOfile = open('range', 'w')
    orbitOfile = open('orbit', 'w')
    sysinfoOfile = open('sysinfo', 'w')

    Header =struct.unpack('H',readFirstHeader)[0]
    #print 'headerID=', hex(Header)
    if Header != 31103:
        print(' \n error with file, check to see that file is a waves file')
        break

    firstLeaderHeader=259
    wavePingHeader=515
    lastLeaderHeader=771
    HPRheader=1027

    #read the checksum offset
    csOffset=ifile.read(2)
    csOffset=struct.unpack('H',csOffset)[0]

    #read 4 bytes before the leader ID
    spare=ifile.read(1)

    readDatatypes=ifile.read(1)
    datatypes=struct.unpack('B',readDatatypes)[0]

    readOffset=ifile.read(2)
    offset=struct.unpack('H',readOffset)[0]

    #this is added in to account for extra bytes in the header
    dataOffsetDif=offset-8
    ifile.read(dataOffsetDif)
    

    #THIS SECTION IS FOR THE FIRST DATA TYPE
    
    
    readLeaderHeader=ifile.read(2)
    leaderHeader =struct.unpack('H',readLeaderHeader)[0]
    #print hex(leaderHeader)
    
    if leaderHeader != firstLeaderHeader:
        print('error with first leader header')
        print hex(leaderHeader)
        break     
    #read and write sysinfo

    #firmware version
    readFwVersion=ifile.read(2)
    fwVersion=struct.unpack('H',readFwVersion)[0]
    sysinfoOfile.write(str(fwVersion))
    sysinfoOfile.write('\n')

    #Configuration (not sure what to do with this)
    ifile.read(2)
        
    # number of depth cells
    readNbins=ifile.read(1)
    Nbins= struct.unpack('B',readNbins)[0]
    sysinfoOfile.write(str(Nbins))
    sysinfoOfile.write('\n')

    #Samples per wave burst
    readNumSamples=ifile.read(2)
    numSamples=struct.unpack('H',readNumSamples)[0]
    sysinfoOfile.write(str(numSamples))
    sysinfoOfile.write('\n')

    # Size of depth cell
    readBinLength=ifile.read(2)
    binLength= struct.unpack('H',readBinLength)[0]
    #change size to meters
    binLength=binLength/100.00
    sysinfoOfile.write(str(binLength))
    sysinfoOfile.write('\n')

    #Time between wave samples
    readTBS=ifile.read(2)
    TBS=struct.unpack('H',readTBS)[0]
    sysinfoOfile.write(str(TBS))
    sysinfoOfile.write('\n')

    #Time between wave bursts
    readTBB=ifile.read(2)
    TBB=struct.unpack('H',readTBB)[0]
    sysinfoOfile.write(str(TBB))
    sysinfoOfile.write('\n')
        
    # distance to first bin
    readDist1=ifile.read(2)
    dist1= struct.unpack('H',readDist1)[0]
    #change distance to meters
    dist1=dist1/100.00
    sysinfoOfile.write(str(dist1))
    sysinfoOfile.write('\n')

        
    # number of bins output
    readBinsOut=ifile.read(1)
    binsOut= struct.unpack('B',readBinsOut)[0]
    sysinfoOfile.write(str(binsOut))
    sysinfoOfile.write('\n')

    #reserved data (?)
    ifile.read(2)

    #skip the directional bitmap
    ifile.read(16)
        
    #read and write the bitmap for what bins are output
    readByteList=ifile.read(16)
    byteList=struct.unpack('16B',readByteList)
        
    def int2msbits(n, count=8):
        """returns the most-significant binary of integer n, using count number of digits"""
        return "".join([str((n >> y) & 1) for y in range(0, count)])

    bytemap = ''
    bitmap = ''
    for ind in range(len(byteList)):
        byte = struct.pack('B', byteList[ind])
        bytemap += byte
        bits = int2msbits(byteList[ind])
        bitmap += bits

    # find indexes from bit map
    # (1) one way to get indices of 1's
    ind = [];
    for i in range(len(bitmap)):
        if bitmap[i]=='1': ind.append(i)

    # transducer height above the bottom
    th = thab

    # heights above transducer
    hat = [dist1+(binLength*elem) for elem in ind]

    # heights above bottom
    hab = [th+elem for elem in hat]

    #if the bins output = 3 or 4 need to include NaNs

    if binsOut == 4:
        sysinfoOfile.write('NaN')
        sysinfoOfile.write('\n')
    elif binsOut == 3:
        sysinfoOfile.write('NaN')
        sysinfoOfile.write('\n')        
        sysinfoOfile.write('NaN')
        sysinfoOfile.write('\n')
        
    for value in hab:
        sysinfoOfile.write(str(value))
        sysinfoOfile.write('\n')


    #Read start time and write to file
    readCentury=ifile.read(1)
    century =struct.unpack('B',readCentury)[0]

    readYear=ifile.read(1)
    year =struct.unpack('B',readYear)[0]

    readMonth=ifile.read(1)
    month=struct.unpack('B',readMonth)[0]

    readDay=ifile.read(1)
    day=struct.unpack('B',readDay)[0]

    readHour=ifile.read(1)
    hour=struct.unpack('B',readHour)[0]

    readMinute=ifile.read(1)
    minute=struct.unpack('B',readMinute)[0]

    readSecond=ifile.read(1)
    second=struct.unpack('B',readSecond)[0]

    readHundSec=ifile.read(1)
    hundSec=struct.unpack('B',readHundSec)[0]

    startTime=str(century)+str('%02d' % year)+str('%02d' % month)+str('%02d' % day)+str('%02d' % hour)+str('%02d' % minute)+str('%02d' % second)+str('%02d' % hundSec)

    #create a time stamp for the file names that only includes from year up to minutes
    FileStamp=str('%02d' % year)+str('%02d' % month)+str('%02d' % day)+str('%02d' % hour)+str('%02d' % minute)

    sysinfoOfile.write(startTime)
    sysinfoOfile.write('\n')

    #Burst number (how to unpack this? dont think we need it)
    readBurstNum=ifile.read(4)

    #Serial number (how to unpack this? dont think we need it)
    readSerial=ifile.read(8)

    #Temperature (dont need this, use the avg temp instead)
    readTemp=ifile.read(2)

    
    #Reserved space (?)
    ifile.read(4)

    cs=ifile.read(2)
    
    #should be the packets ID header again
    test=ifile.read(2) 
    test=struct.unpack('H',test)[0]
    #print hex(test);



    # Done with the First Leader Type>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>


    #Wave Ping Type
    
    leaderHeader = wavePingHeader
    while leaderHeader == wavePingHeader:  

    
        #read the checksum offset
        csOffset=ifile.read(2)
        if len(csOffset) == 0:
            print('\n End of file')
            break
        csOffset=struct.unpack('H',csOffset)[0]
        #print csOffset
        

        #read the spare and datatype
        spare=ifile.read(1)

        datatypes=ifile.read(1)

        #read the offset until the next datatype
        offset=ifile.read(2)
        offset=struct.unpack('H',offset)[0]

        #this is added in to account for extra bytes in the header
        dataOffsetDif=offset-8
        ifile.read(dataOffsetDif)

    
        readLeaderHeader=ifile.read(2)
        leaderHeader =struct.unpack('H',readLeaderHeader)[0]

        # if the leader header isnt the wave ping header move on to the next leader
        if leaderHeader != wavePingHeader:
            break

        #read through sample number and time since burst (don't need this)
        ifile.read(6)

        #read and write the pressure data
        pressure=ifile.read(4)
        pressure=struct.unpack('I',pressure)[0]

        pressureOfile.write(str(pressure))
        pressureOfile.write('\n')

        #read and write the range data

        range1=ifile.read(4)
        range1=struct.unpack('I',range1)[0]
        rangeOfile.write(str(range1))
        rangeOfile.write('  ')

        range2=ifile.read(4)
        range2=struct.unpack('I',range2)[0]
        rangeOfile.write(str(range2))
        rangeOfile.write('  ')

        range3=ifile.read(4)
        range3=struct.unpack('I',range3)[0]
        rangeOfile.write(str(range3))
        rangeOfile.write('  ')

        range4=ifile.read(4)
        range4=struct.unpack('I',range4)[0]
        rangeOfile.write(str(range4))
        rangeOfile.write('\n')


        #read and write the orbital velocity data

        bytes=binsOut*8

        while bytes > 0:
            orbit=ifile.read(2)
            orbit=struct.unpack('h',orbit)[0]
            orbitOfile.write('%06d' % orbit)
            orbitOfile.write('  ')
            bytes=bytes-2

        orbitOfile.write('\n')


    # Done with the Wave Ping Type>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>

        # if the ADCP was set up to take heading, pitch and roll measurements
        # for every ensemble needs to include this

        ID=ifile.read(2)
        ID=struct.unpack('H',ID)[0]
        #print hex(ID)

        if ID==HPRheader:
            heading=ifile.read(2)
            pitch=ifile.read(2)
            roll=ifile.read(2)
            spare=ifile.read(2) #not sure what this is
            cs=ifile.read(2)
        else:
            ifile.seek(-2, 1)
            spare=ifile.read(2)
            cs=ifile.read(2)
        
        
        ID=ifile.read(2)
        if len(ID) == 0:
            break
        ID=struct.unpack('H',ID)[0]
    
        

    
    #readLeaderHeader=ifile.read(2)
    #leaderHeader =struct.unpack('H',readLeaderHeader)[0]
    
    if leaderHeader != lastLeaderHeader:
        break

    #average depth
    readAvgDepth= ifile.read(2)
    avgDepth= struct.unpack('H', readAvgDepth)[0]
    sysinfoOfile.write(str(avgDepth))
    sysinfoOfile.write('\n')

    #Average speed of sound in m/s
    readAvgSound=ifile.read(2)
    avgSound=struct.unpack('H',readAvgSound)[0]
    sysinfoOfile.write(str(avgSound))
    sysinfoOfile.write('\n')

    #Average Temp in .01 deg C
    readAvgTemp=ifile.read(2)
    avgTemp=struct.unpack('H',readAvgTemp)[0]
    sysinfoOfile.write(str(avgTemp))
    sysinfoOfile.write('\n')
        
    #heading
    readHeading= ifile.read(2)
    heading= struct.unpack('H', readHeading)[0]
    sysinfoOfile.write(str(heading))
    sysinfoOfile.write('\n')

    #skip the std dev heading
    ifile.read(2)

    #pitch
    readPitch= ifile.read(2)
    pitch= struct.unpack('h', readPitch)[0]
    sysinfoOfile.write(str(pitch))
    sysinfoOfile.write('\n')   

    #skip the std dev pitch
    stdPitch=ifile.read(2)

    #Roll
    readRoll= ifile.read(2)
    roll= struct.unpack('h', readRoll)[0]
    sysinfoOfile.write(str(roll))
    sysinfoOfile.write('\n')

    #skip std dev roll and last 4 bytes (what are they? maybe a checksum?)
    stdRoll=ifile.read(2)

    cs=ifile.read(2)

    spare=ifile.read(2)
    
    


    pressureOfile.close()
    rangeOfile.close()
    orbitOfile.close()
    sysinfoOfile.close()


    shutil.copy('pressure','pressure_'+FileStamp+'.txt')
    shutil.copy('orbit','orbit_'+FileStamp+'.txt')
    shutil.copy('range','range_'+FileStamp+'.txt')       
    shutil.copy('sysinfo','sysinfo_'+FileStamp+'.txt')

    os.remove('pressure')
    os.remove('orbit')
    os.remove('range')
    os.remove('sysinfo')
    

ifile.close()
print 'processing complete'