Changeset 19432


Ignore:
Timestamp:
07/10/15 06:35:43 (10 years ago)
Author:
cborstad
Message:

CHG: implemented basic framework for quiver plot in python

Location:
issm/trunk-jpl/src/m/plot
Files:
3 edited

Legend:

Unmodified
Added
Removed

  • issm/trunk-jpl/src/m/plot/plot_overlay.py

    r19425 r19432  
    3636        if not is2d:
    3737                raise StandardError('overlay plot not supported for 3D meshes, project on a 2D layer first')
    38         if datatype==3:
    39                 raise StandardError('overlay not yet supported for quiver plots')       
    4038
    4139        if not options.exist('geotiff_name'):

  • issm/trunk-jpl/src/m/plot/plot_unit.py

    r17953 r19432  
    11from cmaptools import truncate_colormap
    22try:
    3         import pylab as p
    4         import matplotlib as mpl
    5         import matplotlib.pyplot as plt
     3    import pylab as p
     4    import matplotlib as mpl
     5    import matplotlib.pyplot as plt
     6    import numpy as npy
    67except ImportError:
    7         print "could not import pylab, matplotlib has not been installed, no plotting capabilities enabled"
     8    print "could not import pylab, matplotlib has not been installed, no plotting capabilities enabled"
    89
    910def plot_unit(x,y,z,elements,data,is2d,isplanet,datatype,options,ax):
    10         """
    11         PLOT_UNIT - unit plot, display data
    12        
    13                 Usage:
    14                         plot_unit(x,y,z,elements,data,is2d,isplanet,datatype,options)
    15        
    16                 See also: PLOTMODEL, PLOT_MANAGER
    17         """
    18        
    19         #edgecolor
    20         edgecolor=options.getfieldvalue('edgecolor','None')
    21        
    22         #number of colorlevels for plots
    23         colorlevels=options.getfieldvalue('colorlevels',128)
    24        
    25         alpha=options.getfieldvalue('alpha',1)
    26        
    27         #colormap
    28         # default sequential colormap
    29         defaultmap=truncate_colormap(mpl.cm.gnuplot2,0.1,0.9,128)
    30         cmap=options.getfieldvalue('colormap',defaultmap)
    31         if options.exist('cmap_set_over'):
    32            over=options.getfieldvalue('cmap_set_over','0.5')
    33            cmap.set_over(over)
    34         if options.exist('cmap_set_under'):
    35            under=options.getfieldvalue('cmap_set_under','0.5')
    36            cmap.set_under(under)
    37        
    38         #normalize colormap if clim/caxis specified
    39         if options.exist('clim'):
    40            lims=options.getfieldvalue('clim',[min(data),max(data)])
    41         elif options.exist('caxis'):
    42            lims=options.getfieldvalue('caxis',[min(data),max(data)])
    43         else:
    44                 if min(data)==max(data):
    45                         lims=[min(data)-0.5,max(data)+0.5]
    46                 else:
    47                         lims=[min(data),max(data)]
    48         norm = mpl.colors.Normalize(vmin=lims[0], vmax=lims[1])
    49         if datatype==1:
    50            #element plot
    51                 if is2d:
    52                         tri=ax.tripcolor(x,y,elements,data,colorlevels,cmap=cmap,norm=norm,alpha=alpha,edgecolors=edgecolor)
    53                 else:
    54                         raise ValueError('plot_unit error: 3D element plot not supported yet')
    55                 return
    56        
    57         elif datatype==2:
    58                 #node plot
    59                 if is2d:
    60                         tri=ax.tricontourf(x,y,elements,data,colorlevels,cmap=cmap,norm=norm,alpha=alpha,extend='both')
    61                         if edgecolor != 'None':
    62                                 ax.triplot(x,y,elements,color=edgecolor)
    63                 else:
    64                         raise ValueError('plot_unit error: 3D node plot not supported yet')
    65                 return
    66        
    67         elif datatype==3:
    68                 print 'plot_unit message: quiver plot not implemented yet'
    69                 return
    70        
    71         elif datatype==4:
    72                 #P1 patch plot
    73                 print 'plot_unit message: P1 patch plot not implemented yet'
    74                 return
    75        
    76         elif datatype==5:
    77                 print 'plot_unit message: P0 patch plot not implemented yet'
    78                 return
    79        
    80         else:
    81                 raise ValueError('datatype=%d not supported' % datatype)
    82 
     11    """
     12    PLOT_UNIT - unit plot, display data
     13   
     14        Usage:
     15                plot_unit(x,y,z,elements,data,is2d,isplanet,datatype,options)
     16   
     17        See also: PLOTMODEL, PLOT_MANAGER
     18    """
     19   
     20    #edgecolor
     21    edgecolor=options.getfieldvalue('edgecolor','None')
     22   
     23    #number of colorlevels for plots
     24    colorlevels=options.getfieldvalue('colorlevels',128)
     25   
     26    alpha=options.getfieldvalue('alpha',1)
     27   
     28    #colormap
     29    # default sequential colormap
     30    defaultmap=truncate_colormap(mpl.cm.gnuplot2,0.1,0.9,128)
     31    cmap=options.getfieldvalue('colormap',defaultmap)
     32    if options.exist('cmap_set_over'):
     33        over=options.getfieldvalue('cmap_set_over','0.5')
     34        cmap.set_over(over)
     35    if options.exist('cmap_set_under'):
     36        under=options.getfieldvalue('cmap_set_under','0.5')
     37        cmap.set_under(under)
     38   
     39    #normalize colormap if clim/caxis specified
     40    if options.exist('clim'):
     41        lims=options.getfieldvalue('clim',[npy.amin(data),npy.amax(data)])
     42    elif options.exist('caxis'):
     43        lims=options.getfieldvalue('caxis',[npy.amin(data),npy.amax(data)])
     44    else:
     45        if npy.amin(data)==npy.amax(data):
     46            lims=[npy.amin(data)-0.5,npy.amax(data)+0.5]
     47        else:
     48            lims=[npy.amin(data),npy.amax(data)]
     49    norm = mpl.colors.Normalize(vmin=lims[0], vmax=lims[1])
     50    if datatype==1:
     51       #element plot
     52        if is2d:
     53            tri=ax.tripcolor(x,y,elements,data,colorlevels,cmap=cmap,norm=norm,alpha=alpha,edgecolors=edgecolor)
     54        else:
     55            raise ValueError('plot_unit error: 3D element plot not supported yet')
     56        return
     57   
     58    elif datatype==2:
     59        #node plot
     60        if is2d:
     61            tri=ax.tricontourf(x,y,elements,data,colorlevels,cmap=cmap,norm=norm,alpha=alpha,extend='both')
     62            if edgecolor != 'None':
     63                ax.triplot(x,y,elements,color=edgecolor)
     64        else:
     65            raise ValueError('plot_unit error: 3D node plot not supported yet')
     66        return
     67   
     68    elif datatype==3:
     69        vx=data[:,0]
     70        vy=data[:,1]
     71        #TODO write plot_quiver.py to handle this here
     72        color=npy.sqrt(vx**2+vy**2)
     73        scale=options.getfieldvalue('scale',1000)
     74        width=options.getfieldvalue('width',0.005*(npy.amax(x)-npy.amin(y)))
     75        headwidth=options.getfieldvalue('headwidth',3)
     76        headlength=options.getfieldvalue('headlength',5)
     77        Q=ax.quiver(x,y,vx,vy,color,cmap=cmap,norm=norm,scale=scale,
     78                width=width,headwidth=headwidth,headlength=headlength)
     79        return
     80   
     81    elif datatype==4:
     82        #P1 patch plot
     83        print 'plot_unit message: P1 patch plot not implemented yet'
     84        return
     85   
     86    elif datatype==5:
     87        print 'plot_unit message: P0 patch plot not implemented yet'
     88        return
     89   
     90    else:
     91        raise ValueError('datatype=%d not supported' % datatype)
     92   

  • issm/trunk-jpl/src/m/plot/processdata.py

    r19430 r19432  
    33
    44def processdata(md,data,options):
    5         """
    6         PROCESSDATA - process data to be plotted
    7 
    8                 datatype = 1 -> elements
    9                 datatype = 2 -> nodes
    10                 datatype = 3 -> node quivers
    11                 datatype = 4 -> patch
    12 
    13                 Usage:
    14                         data,datatype=processdata(md,data,options);
    15 
    16                 See also: PLOTMODEL, PROCESSMESH
    17         """
    18 
    19         #check format
    20         if (len(data)==0 or (len(data)==1 and not isinstance(data,dict) and isnan(data).all())):
    21                 raise ValueError("processdata error message: 'data' provided is empty")
    22 
    23         #needed later on
    24         if 'numberofvertices2d' in dir(md.mesh):
    25                 numberofvertices2d=md.mesh.numberofvertices2d
    26                 numberofelements2d=md.mesh.numberofelements2d
    27         else:
    28                 numberofvertices2d=npy.nan
    29                 numberofelements2d=npy.nan
    30 
    31         procdata=npy.copy(data)
    32 
    33         #process patch
    34 
    35         #initialize datatype
    36         datatype=0
    37 
    38         #get datasize
    39         if npy.ndim(procdata)==1:
    40                 datasize=npy.array([len(procdata),1])
    41         else:
    42                 datasize=npy.shape(procdata)
    43        
    44         #process NaN's if any
    45         nanfill=options.getfieldvalue('nan',-9999)
    46         if npy.any(npy.isnan(procdata)):
    47                 lb=npy.min(data[~npy.isnan(data)])
    48                 ub=npy.max(data[~npy.isnan(data)])
    49                 if lb==ub:
    50                         lb=lb-0.5
    51                         ub=ub+0.5
    52                         nanfill=lb-1
    53                 procdata[npy.isnan(procdata)]=nanfill
    54                 options.addfielddefault('clim',[lb,ub])
    55                 options.addfielddefault('cmap_set_under','1')
    56                 print "WARNING: nan's treated as", nanfill, "by default.  Change using pairoption 'nan',nan_fill_value in plotmodel call"
    57 
    58         #non-patch processing
    59 
    60         #element data
    61         if datasize[0]==md.mesh.numberofelements and datasize[1]==1:
    62                
    63                 #initialize datatype if non patch
    64                 if datatype!=4 and datatype!=5:
    65                         datatype=1
    66 
    67                 #mask?
    68 
    69                 #log?
    70 
    71         #node data
    72         if datasize[0]==md.mesh.numberofvertices and datasize[1]==1:
    73                 datatype=2
    74 
    75         #spc time series?
    76         if datasize[0]==md.mesh.numberofvertices+1:
    77                 datatype=2
    78                 spccol=options.getfieldvalue('spccol',0)
    79                 print 'multiple-column spc field; specify column to plot using option "spccol"'
    80                 print 'column ', spccol, ' plotted for time: ', procdata[-1,spccol]
    81                 procdata=procdata[0:-1,spccol]
    82 
    83                 #mask?
    84 
    85                 #log?
    86 
    87         #layer projection?
    88 
    89         #control arrow density if quiver plot
    90 
    91         #convert rank-2 array to rank-1
    92         if npy.ndim(procdata)==2 and npy.shape(procdata)[1]==1:
    93                 procdata=procdata.reshape(-1,)
    94 
    95         #if datatype is still zero, error out
    96         if datatype==0:
    97                 raise ValueError("processdata error: data provided not recognized or not supported")
    98         else:
    99                 return procdata, datatype
     5    """
     6    PROCESSDATA - process data to be plotted
     7   
     8        datatype = 1 -> elements
     9        datatype = 2 -> nodes
     10        datatype = 3 -> node quivers
     11        datatype = 4 -> patch
     12   
     13        Usage:
     14                data,datatype=processdata(md,data,options);
     15   
     16        See also: PLOTMODEL, PROCESSMESH
     17    """
     18   
     19    #check format
     20    if (len(data)==0 or (len(data)==1 and not isinstance(data,dict) and isnan(data).all())):
     21        raise ValueError("processdata error message: 'data' provided is empty")
     22   
     23    #needed later on
     24    if 'numberofvertices2d' in dir(md.mesh):
     25        numberofvertices2d=md.mesh.numberofvertices2d
     26        numberofelements2d=md.mesh.numberofelements2d
     27    else:
     28        numberofvertices2d=npy.nan
     29        numberofelements2d=npy.nan
     30   
     31    procdata=npy.copy(data)
     32   
     33    #process patch
     34   
     35    #initialize datatype
     36    datatype=0
     37   
     38    #get datasize
     39    if npy.ndim(procdata)==1:
     40        datasize=npy.array([len(procdata),1])
     41    else:
     42        datasize=npy.shape(procdata)
     43        if len(datasize)>2:
     44            raise ValueError('data passed to plotmodel has more than 2 dimensions; check that column vectors are rank-1')
     45   
     46    #process NaN's if any
     47    nanfill=options.getfieldvalue('nan',-9999)
     48    if npy.any(npy.isnan(procdata)):
     49        lb=npy.min(data[~npy.isnan(data)])
     50        ub=npy.max(data[~npy.isnan(data)])
     51        if lb==ub:
     52            lb=lb-0.5
     53            ub=ub+0.5
     54            nanfill=lb-1
     55        procdata[npy.isnan(procdata)]=nanfill
     56        options.addfielddefault('clim',[lb,ub])
     57        options.addfielddefault('cmap_set_under','1')
     58        print "WARNING: nan's treated as", nanfill, "by default.  Change using pairoption 'nan',nan_fill_value in plotmodel call"
     59   
     60    #quiver plot
     61    if datasize[1]>1:
     62        if datasize[0]==md.mesh.numberofvertices and datasize[1]==2:
     63            datatype=3
     64        else:
     65            raise ValueError('plotmodel error message: data should have two columns of length md.mesh.numberofvertices for a quiver plot')
     66   
     67    #non-patch processing
     68   
     69    #element data
     70    if datasize[0]==md.mesh.numberofelements and datasize[1]==1:
     71       
     72        #initialize datatype if non patch
     73        if datatype!=4 and datatype!=5:
     74            datatype=1
     75   
     76        #mask?
     77   
     78        #log?
     79   
     80    #node data
     81    if datasize[0]==md.mesh.numberofvertices and datasize[1]==1:
     82        datatype=2
     83   
     84    #spc time series?
     85    if datasize[0]==md.mesh.numberofvertices+1:
     86        datatype=2
     87        spccol=options.getfieldvalue('spccol',0)
     88        print 'multiple-column spc field; specify column to plot using option "spccol"'
     89        print 'column ', spccol, ' plotted for time: ', procdata[-1,spccol]
     90        procdata=procdata[0:-1,spccol]
     91   
     92        #mask?
     93   
     94        #log?
     95   
     96    #layer projection?
     97   
     98    #control arrow density if quiver plot
     99   
     100    #convert rank-2 array to rank-1
     101    if npy.ndim(procdata)==2 and npy.shape(procdata)[1]==1:
     102        procdata=procdata.reshape(-1,)
     103   
     104    #if datatype is still zero, error out
     105    if datatype==0:
     106        raise ValueError("processdata error: data provided not recognized or not supported")
     107    else:
     108        return procdata, datatype
Note: See TracChangeset for help on using the changeset viewer.