root/sodar/trunk/plotSodar.py

#!/usr/bin/python
"""
Module to process sodar data.
"""

import optparse
from sodar.utils import openAnything, findMissing 
from sodar import arrayData as a
import numpy as n
import pylab as p
import os

def plotSingle(filein, pathout):
    """
    Plot singe sodar raw data file.
    
    plotSingle(filein, pathout) -> rc
    
    filein - path to raw data file.
    fileout - path to plot images directory.
    rc - 0 if successful
    """
    
    if not os.path.exists(pathout):
        os.makedirs(pathout, mode=0775)
    
    try:
        rawDataHandle = openAnything.openAnything(filein)
        rawDataString = rawDataHandle.read()
        rawDataHandle.close()
    except:
        raise IOError("Failure to read raw data.")
    
    rawDataObject = a.rawData.RawData(rawDataString)
    formattedDataObject = a.formattedData.FormattedData(rawDataObject)
    arrayDataObject = a.ArrayData(formattedDataObject)
    
    beginStamp = arrayDataObject.beginStamp
    endStamp = arrayDataObject.endStamp
    numIntervals = arrayDataObject.numIntervals
    timeInterval = arrayDataObject.timeInterval
    timeInterval = (((float(timeInterval.days) * 24.0) +
                     (timeInterval.seconds / 3600.0)) * 60.0)
    minAltitude = arrayDataObject.minAltitude
    altInterval = arrayDataObject.altInterval
    
    def makePcolor(vector, title, fileout, cbLabel):
        fig = p.figure(1)
        axe = fig.add_subplot(1, 1, 1)
        pc = axe.pcolor(vector)
        
        axe.set_xlabel('Time (min)')
        axe.set_ylabel('Altitude (m)')
        axe.set_xbound(upper=numIntervals)
        
        xticks = axe.get_xticks()
        xticklabels = [str(x * timeInterval)
                       for x in xticks]
        axe.set_xticklabels(xticklabels)
        
        yticks = axe.get_yticks()
        yticklabels = [str(y * altInterval + minAltitude)
                       for y in yticks]
        axe.set_yticklabels(yticklabels)
        
        axe.set_title(title)
        cb = p.colorbar(pc)
        cb.set_label(cbLabel)
        
        fig.savefig(os.path.join(pathout, fileout))
        fig.clear()
        
        return 0
    
    uComponents = arrayDataObject.uComponents
    maskedUComponents = n.ma.masked_where(n.isnan(uComponents),
                                          uComponents)
    makePcolor(maskedUComponents.T,
               'U Component of Wind Veloctiy\nFrom %s To %s' %
               (beginStamp, endStamp),
               'uComponents.png',
               'Speed (cm/sec)')
    
    vComponents = arrayDataObject.vComponents
    maskedVComponents = n.ma.masked_where(n.isnan(vComponents),
                                          vComponents)
    makePcolor(maskedVComponents.T,
               'V Component of Wind Veloctiy\nFrom %s To %s' %
               (beginStamp, endStamp),
               'vComponents.png',
               'Speed (cm/sec)')
    
    wComponents = arrayDataObject.wComponents
    maskedWComponents = n.ma.masked_where(n.isnan(wComponents),
                                          wComponents)
    makePcolor(maskedWComponents.T,
               'W Component of Wind Veloctiy\nFrom %s To %s' %
               (beginStamp, endStamp),
               'wComponents.png',
               'Speed (cm/sec)')
    
    echoStrengths = arrayDataObject.echoStrengths
    maskedEchoStrenths = n.ma.masked_where(n.isnan(echoStrengths),
                                                   echoStrengths)
    makePcolor(maskedEchoStrenths.T,
               'Echo Strength\nFrom %s To %s' %
               (beginStamp, endStamp),
               'echoStrengths.png',
               'Strength (no units)')
    
    fig = p.figure(1)
    axe = fig.add_subplot(1, 1, 1)
    qv = axe.quiver(maskedUComponents,
                    maskedVComponents,
                    maskedWComponents)
    axe.set_xlabel('Time (min)')
    axe.set_ylabel('Altitude (m)')
    axe.set_xbound(upper=numIntervals)
    xticks = axe.get_xticks()
    xticklabels = [str(x * timeInterval)
                   for x in xticks]
    axe.set_xticklabels(xticklabels)
    yticks = axe.get_yticks()
    yticklabels = [str(y * altInterval + minAltitude)
                   for y in yticks]
    axe.set_yticklabels(yticklabels)
    axe.set_title('Wind Velocty\nFrom %s To %s' %
                  (beginStamp, endStamp))
    cb = p.colorbar(qv)
    cb.set_label('W Component Speed (cm/sec)')
    fig.savefig(os.path.join(pathout, 'quiver.png'))
    fig.clear()
    
    return 0

def plotAll(source, destination):
    """
    Plot all sodar raw data files.
    
    plotAll(source, destination) -> rc
    
    source - path to raw data directory in NCCOOS format.
    detination - path to plot images directory in NCCOOS format.
    rc - 0 if successful
    """
    
    if source.endswith(os.path.sep):
        source = source.rstrip(os.path.sep)
    
    if destination.endswith(os.path.sep):
        destination = destination.rstrip(os.path.sep)
    
    walkList = findMissing.findMissing(source, destination)
    
    for filein, pathout in walkList:
        plotSingle(filein, pathout)
    
    return 0

def _main():
    """bin/python %prog [options] /path/to/raw/files/ /path/to/plots/dir/"""
    
    __description__ = 'Plot all sodar raw data files.'
    
    parser = optparse.OptionParser(usage=_main.__doc__,
                                   version='1.0',
                                   description=__description__)
    (values, args) = parser.parse_args()
    (source, destination) = tuple(args)
    plotAll(source, destination)
    
    return 0

if __name__ == "__main__":
    _main()