Changeset 23920


Ignore:
Timestamp:
05/14/19 13:16:58 (6 years ago)
Author:
kruegern
Message:

NEW: added colormap support to python plotmodel(), all matplotlib colormaps, Rignot, and demmap(); default is viridis

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

Legend:

Unmodified
Added
Removed

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

    r23716 r23920  
    11import numpy as  np
    2 from cmaptools import truncate_colormap
     2from cmaptools import getcolormap
    33from plot_contour import plot_contour
    44from plot_streamlines import plot_streamlines
    … …  
    66
    77try:
    8         from matplotlib.ticker import MaxNLocator
    9         from mpl_toolkits.axes_grid1 import make_axes_locatable
    10         from mpl_toolkits.mplot3d import Axes3D
    11         import matplotlib as mpl
    12         import matplotlib.pyplot as plt
     8    from matplotlib.ticker import MaxNLocator
     9    from mpl_toolkits.axes_grid1 import make_axes_locatable
     10    from mpl_toolkits.mplot3d import Axes3D
     11    import matplotlib as mpl
     12    import matplotlib.pyplot as plt
    1313except ImportError:
    14         print("could not import pylab, matplotlib has not been installed, no plotting capabilities enabled")
     14    print("could not import pylab, matplotlib has not been installed, no plotting capabilities enabled")
    1515
    1616def applyoptions(md,data,options,fig,axgrid,gridindex):
    17         '''
    18         APPLYOPTIONS - apply options to current plot
    19 
    20         'plotobj' is the object returned by the specific plot call used to
    21         render the data.  This object is used for adding a colorbar.
    22 
    23                 Usage:
    24                         applyoptions(md,data,options)
    25 
    26                 See also: PLOTMODEL, PARSE_OPTIONS
    27         '''
    28 
    29         # get handle to current figure and axes instance
    30         #fig = p.gcf()
    31         ax=     axgrid[gridindex]
    32 
    33         # {{{ font
    34         fontsize=options.getfieldvalue('fontsize',8)
    35         fontweight=options.getfieldvalue('fontweight','normal')
    36         fontfamily=options.getfieldvalue('fontfamily','sans-serif')
    37         font={'fontsize'                :fontsize,
    38                                 'fontweight'    :fontweight,
    39                                 'family'                        :fontfamily}
    40         # }}}
    41         # {{{ title
    42         if options.exist('title'):
    43                 title=options.getfieldvalue('title')
    44                 if options.exist('titlefontsize'):
    45                         titlefontsize=options.getfieldvalue('titlefontsize')
    46                 else:
    47                         titlefontsize=fontsize
    48                 if options.exist('titlefontweight'):
    49                         titlefontweight=options.getfieldvalue('titlefontweight')
    50                 else:
    51                         titlefontweight=fontweight
    52                 #title font
    53                 titlefont=font.copy()
    54                 titlefont['size']=titlefontsize
    55                 titlefont['weight']=titlefontweight
    56                 ax.set_title(title,**titlefont)
    57         # }}}
    58         # {{{ xlabel, ylabel, zlabel
    59         if options.exist('labelfontsize'):
    60                 labelfontsize=options.getfieldvalue('labelfontsize')
    61         else:
    62                 labelfontsize=fontsize
    63         if options.exist('labelfontweight'):
    64                 labelfontweight=options.getfieldvalue('labelfontweight')
    65         else:
    66                 labelfontweight=fontweight
    67 
    68         #font dict for labels
    69         labelfont=font.copy()
    70         labelfont['fontsize']=labelfontsize
    71         labelfont['fontweight']=labelfontweight
    72 
    73         if options.exist('xlabel'):
    74                 ax.set_xlabel(options.getfieldvalue('xlabel'),**labelfont)
    75         if options.exist('ylabel'):
    76                 ax.set_ylabel(options.getfieldvalue('ylabel'),**labelfont)
    77         if options.exist('zlabel'):
    78                 ax.set_zlabel(options.getfieldvalue('zlabel'),**labelfont)
    79         # }}}
    80         # {{{ xticks, yticks, zticks (tick locations)
    81         if options.exist('xticks'):
    82                 if options.exist('xticklabels'):
    83                         xticklabels=options.getfieldvalue('xticklabels')
    84                         ax.set_xticks(options.getfieldvalue('xticks'),xticklabels)
    85                 else:
    86                         ax.set_xticks(options.getfieldvalue('xticks'))
    87         if options.exist('yticks'):
    88                 if options.exist('yticklabels'):
    89                         yticklabels=options.getfieldvalue('yticklabels')
    90                         ax.set_yticks(options.getfieldvalue('yticks'),yticklabels)
    91                 else:
    92                         ax.set_yticks(options.getfieldvalue('yticks'))
    93         if options.exist('zticks'):
    94                 if options.exist('zticklabels'):
    95                         zticklabels=options.getfieldvalue('zticklabels')
    96                         ax.set_zticks(options.getfieldvalue('zticks'),zticklabels)
    97                 else:
    98                         ax.set_zticks(options.getfieldvalue('zticks'))
    99         # }}}
    100         # {{{ xticklabels,yticklabels,zticklabels
    101         if options.getfieldvalue('ticklabels','off')=='off' or options.getfieldvalue('ticklabels',0)==0:
    102                 options.addfielddefault('xticklabels',[])
    103                 options.addfielddefault('yticklabels',[])
    104                 # TODO check if ax has a z-axis (e.g. is 3D)
    105         if options.exist('xticklabels'):
    106                 xticklabels=options.getfieldvalue('xticklabels')
    107                 ax.set_xticklabels(xticklabels)
    108         if options.exist('yticklabels'):
    109                 yticklabels=options.getfieldvalue('yticklabels')
    110                 ax.set_yticklabels(yticklabels)
    111         if options.exist('zticklabels'):
    112                 zticklabels=options.getfieldvalue('zticklabels')
    113                 ax.set_zticklabels(zticklabels)
    114         # }}}
    115         # {{{ ticklabel notation
    116         #ax.ticklabel_format(style='sci',scilimits=(0,0))
    117         # }}}
    118         # {{{ ticklabelfontsize
    119         if options.exist('ticklabelfontsize'):
    120                 for label in ax.get_xticklabels() + ax.get_yticklabels():
    121                         label.set_fontsize(options.getfieldvalue('ticklabelfontsize'))
    122                 if int(md.mesh.dimension)==3:
    123                         for label in ax.get_zticklabels():
    124                                 label.set_fontsize(options.getfieldvalue('ticklabelfontsize'))
    125         # }}}
    126         # {{{ view TOFIX
    127         #if int(md.mesh.dimension) == 3 and options.exist('layer'):
    128         #       #options.getfieldvalue('view') ?
    129         #       ax=fig.gca(projection='3d')
    130         #plt.show()
    131         # }}}
    132         # {{{ axis
    133         if options.exist('axis'):
    134                 if options.getfieldvalue('axis',True)=='off':
    135                         ax.ticklabel_format(style='plain')
    136                         p.setp(ax.get_xticklabels(), visible=False)
    137                         p.setp(ax.get_yticklabels(), visible=False)
    138         # }}}
    139         # {{{ box
    140         if options.exist('box'):
    141                 eval(options.getfieldvalue('box'))
    142         # }}}
    143         # {{{ xlim, ylim, zlim
    144         if options.exist('xlim'):
    145                 ax.set_xlim(options.getfieldvalue('xlim'))
    146         if options.exist('ylim'):
    147                 ax.set_ylim(options.getfieldvalue('ylim'))
    148         if options.exist('zlim'):
    149                 ax.set_zlim(options.getfieldvalue('zlim'))
    150         # }}}
    151         # {{{ latlon TODO
    152         # }}}
    153         # {{{ Basinzoom TODO
    154         # }}}
    155         # {{{ ShowBasins TODO
    156         # }}}
    157         # {{{ clim
    158         if options.exist('clim'):
    159                 lims=options.getfieldvalue('clim')
    160                 assert len(lims)==2, 'error, clim should be passed as a list of length 2'
    161         elif options.exist('caxis'):
    162                 lims=options.getfieldvalue('caxis')
    163                 assert len(lims)==2, 'error, caxis should be passed as a list of length 2'
    164                 options.addfielddefault('clim',lims)
    165         else:
    166                 if len(data)>0:
    167                         lims=[data.min(),data.max()]
    168                 else:
    169                         lims=[0,1]
    170         # }}}
    171         # {{{ shading TODO
    172         #if options.exist('shading'):
    173         # }}}
    174         # {{{ grid
    175         if options.exist('grid'):
    176                 if 'on' in options.getfieldvalue('grid','on'):
    177                         ax.grid()
    178         # }}}
    179         # {{{ colormap
    180         if options.exist('colornorm'):
    181                 norm=options.getfieldvalue('colornorm')
    182         if options.exist('colormap'):
    183                 cmap=options.getfieldvalue('colormap')
    184         cbar_extend=0
    185         if options.exist('cmap_set_over'):
    186                 cbar_extend+=1
    187         if options.exist('cmap_set_under'):
    188                 cbar_extend+=2
    189         # }}}
    190         # {{{ contours
    191         if options.exist('contourlevels'):
    192                 plot_contour(md,data,options,ax)
    193         # }}}
    194         # {{{ wrapping TODO
    195         # }}}
    196         # {{{ colorbar
    197         if options.getfieldvalue('colorbar',1)==1:
    198                 if cbar_extend==0:
    199                         extend='neither'
    200                 elif cbar_extend==1:
    201                         extend='max'
    202                 elif cbar_extend==2:
    203                         extend='min'
    204                 elif cbar_extend==3:
    205                         extend='both'
    206                 cb = mpl.colorbar.ColorbarBase(ax.cax,cmap=cmap, norm=norm, extend=extend)
    207                 if options.exist('alpha'):
    208                         cb.set_alpha(options.getfieldvalue('alpha'))
    209                 if options.exist('colorbarnumticks'):
    210                         cb.locator=MaxNLocator(nbins=options.getfieldvalue('colorbarnumticks',5))
    211                 else:
    212                         cb.locator=MaxNLocator(nbins=5) # default 5 ticks
    213                 if options.exist('colorbartickspacing'):
    214                         locs=np.arange(lims[0],lims[1]+1,options.getfieldvalue('colorbartickspacing'))
    215                         cb.set_ticks(locs)
    216                 if options.exist('colorbarlines'):
    217                         locs=np.arange(lims[0],lims[1]+1,options.getfieldvalue('colorbarlines'))
    218                         cb.add_lines(locs,['k' for i in range(len(locs))],np.ones_like(locs))
    219                 if options.exist('colorbarlineatvalue'):
    220                         locs=options.getfieldvalue('colorbarlineatvalue')
    221                         colors=options.getfieldvalue('colorbarlineatvaluecolor',['k' for i in range (len(locs))])
    222                         widths=options.getfieldvalue('colorbarlineatvaluewidth',np.ones_like(locs))
    223                         cb.add_lines(locs,colors,widths)
    224                 if options.exist('colorbartitle'):
    225                         if options.exist('colorbartitlepad'):
    226                                 cb.set_label(options.getfieldvalue('colorbartitle'),
    227                                                                                  labelpad=options.getfieldvalue('colorbartitlepad'),fontsize=fontsize)
    228                         else:
    229                                 cb.set_label(options.getfieldvalue('colorbartitle'),fontsize=fontsize)
    230                 cb.ax.tick_params(labelsize=fontsize)
    231                 cb.solids.set_rasterized(True)
    232                 cb.update_ticks()
    233                 cb.draw_all()
    234                 if options.exist('colorbarfontsize'):
    235                         colorbarfontsize=options.getfieldvalue('colorbarfontsize')
    236                         cb.ax.tick_params(labelsize=colorbarfontsize)
    237                         # cb.set_ticks([0,-10])
    238                         # cb.set_ticklabels([-10,0,10])
    239                 if options.exist('colorbarticks'):
    240                         colorbarticks=options.getfieldvalue('colorbarticks')
    241                         cb.set_ticks(colorbarticks)
    242                 plt.sca(ax) # return to original axes control
    243         # }}}
    244         # {{{ expdisp
    245         if options.exist('expdisp'):
    246                 expdisp(ax,options)
    247         # }}}
    248         # {{{ area TODO
    249         # }}}
    250         # {{{ text
    251         if options.exist('text'):
    252                 text=options.getfieldvalue('text')
    253                 textx=options.getfieldvalue('textx')
    254                 texty=options.getfieldvalue('texty')
    255                 textcolor=options.getfieldvalue('textcolor')
    256                 textweight=options.getfieldvalue('textweight')
    257                 textrotation=options.getfieldvalue('textrotation')
    258                 textfontsize=options.getfieldvalue('textfontsize')
    259                 for label,x,y,size,color,weight,rotation in zip(text,textx,texty,textfontsize,textcolor,textweight,textrotation):
    260                         ax.text(x,y,label,transform=ax.transAxes,fontsize=size,color=color,weight=weight,rotation=rotation)
    261         # }}}
    262         # {{{ north arrow TODO
    263         # }}}
    264         # {{{ scale ruler TODO
    265         # }}}
    266         # {{{ streamlines TOFIX
    267         if options.exist('streamlines'):
    268                 plot_streamlines(md,options,ax)
    269         # }}}
    270         # {{{ axis positions TODO
    271         # }}}
    272         # {{{ figure position TODO
    273         # }}}
    274         # {{{ axes position TODO
    275         # }}}
    276         # {{{ showregion TODO
    277         # }}}
    278         # {{{ flat edges of a partition TODO
    279         # }}}
    280         # {{{ scatter TODO
    281         # }}}
    282         # {{{ backgroundcolor TODO
    283         # }}}
    284         # {{{ figurebackgroundcolor TODO
    285         # }}}
    286         # {{{ lighting TODO
    287         # }}}
    288         # {{{ point cloud TODO
    289         # }}}
    290         # {{{ inset TODO
    291         # }}}
     17    '''
     18    APPLYOPTIONS - apply options to current plot
     19
     20    'plotobj' is the object returned by the specific plot call used to
     21    render the data.  This object is used for adding a colorbar.
     22
     23        Usage:
     24            applyoptions(md,data,options)
     25
     26        See also: PLOTMODEL, PARSE_OPTIONS
     27    '''
     28
     29    # get handle to current figure and axes instance
     30    #fig = p.gcf()
     31    ax=    axgrid[gridindex]
     32
     33    # {{{ font
     34    fontsize=options.getfieldvalue('fontsize',8)
     35    fontweight=options.getfieldvalue('fontweight','normal')
     36    fontfamily=options.getfieldvalue('fontfamily','sans-serif')
     37    font={'fontsize'        :fontsize,
     38                'fontweight'    :fontweight,
     39                'family'            :fontfamily}
     40    # }}}
     41    # {{{ title
     42    if options.exist('title'):
     43        title=options.getfieldvalue('title')
     44        if options.exist('titlefontsize'):
     45            titlefontsize=options.getfieldvalue('titlefontsize')
     46        else:
     47            titlefontsize=fontsize
     48        if options.exist('titlefontweight'):
     49            titlefontweight=options.getfieldvalue('titlefontweight')
     50        else:
     51            titlefontweight=fontweight
     52        #title font
     53        titlefont=font.copy()
     54        titlefont['size']=titlefontsize
     55        titlefont['weight']=titlefontweight
     56        ax.set_title(title,**titlefont)
     57    # }}}
     58    # {{{ xlabel, ylabel, zlabel
     59    if options.exist('labelfontsize'):
     60        labelfontsize=options.getfieldvalue('labelfontsize')
     61    else:
     62        labelfontsize=fontsize
     63    if options.exist('labelfontweight'):
     64        labelfontweight=options.getfieldvalue('labelfontweight')
     65    else:
     66        labelfontweight=fontweight
     67
     68    #font dict for labels
     69    labelfont=font.copy()
     70    labelfont['fontsize']=labelfontsize
     71    labelfont['fontweight']=labelfontweight
     72
     73    if options.exist('xlabel'):
     74        ax.set_xlabel(options.getfieldvalue('xlabel'),**labelfont)
     75    if options.exist('ylabel'):
     76        ax.set_ylabel(options.getfieldvalue('ylabel'),**labelfont)
     77    if options.exist('zlabel'):
     78        ax.set_zlabel(options.getfieldvalue('zlabel'),**labelfont)
     79    # }}}
     80    # {{{ xticks, yticks, zticks (tick locations)
     81    if options.exist('xticks'):
     82        if options.exist('xticklabels'):
     83            xticklabels=options.getfieldvalue('xticklabels')
     84            ax.set_xticks(options.getfieldvalue('xticks'),xticklabels)
     85        else:
     86            ax.set_xticks(options.getfieldvalue('xticks'))
     87    if options.exist('yticks'):
     88        if options.exist('yticklabels'):
     89            yticklabels=options.getfieldvalue('yticklabels')
     90            ax.set_yticks(options.getfieldvalue('yticks'),yticklabels)
     91        else:
     92            ax.set_yticks(options.getfieldvalue('yticks'))
     93    if options.exist('zticks'):
     94        if options.exist('zticklabels'):
     95            zticklabels=options.getfieldvalue('zticklabels')
     96            ax.set_zticks(options.getfieldvalue('zticks'),zticklabels)
     97        else:
     98            ax.set_zticks(options.getfieldvalue('zticks'))
     99    # }}}
     100    # {{{ xticklabels,yticklabels,zticklabels
     101    if options.getfieldvalue('ticklabels','off')=='off' or options.getfieldvalue('ticklabels',0)==0:
     102        options.addfielddefault('xticklabels',[])
     103        options.addfielddefault('yticklabels',[])
     104        # TODO check if ax has a z-axis (e.g. is 3D)
     105    if options.exist('xticklabels'):
     106        xticklabels=options.getfieldvalue('xticklabels')
     107        ax.set_xticklabels(xticklabels)
     108    if options.exist('yticklabels'):
     109        yticklabels=options.getfieldvalue('yticklabels')
     110        ax.set_yticklabels(yticklabels)
     111    if options.exist('zticklabels'):
     112        zticklabels=options.getfieldvalue('zticklabels')
     113        ax.set_zticklabels(zticklabels)
     114    # }}}
     115    # {{{ ticklabel notation
     116    #ax.ticklabel_format(style='sci',scilimits=(0,0))
     117    # }}}
     118    # {{{ ticklabelfontsize
     119    if options.exist('ticklabelfontsize'):
     120        for label in ax.get_xticklabels() + ax.get_yticklabels():
     121            label.set_fontsize(options.getfieldvalue('ticklabelfontsize'))
     122        if int(md.mesh.dimension)==3:
     123            for label in ax.get_zticklabels():
     124                label.set_fontsize(options.getfieldvalue('ticklabelfontsize'))
     125    # }}}
     126    # {{{ view TOFIX
     127    #if int(md.mesh.dimension) == 3 and options.exist('layer'):
     128    #    #options.getfieldvalue('view') ?
     129    #    ax=fig.gca(projection='3d')
     130    #plt.show()
     131    # }}}
     132    # {{{ axis
     133    if options.exist('axis'):
     134        if options.getfieldvalue('axis',True)=='off':
     135            ax.ticklabel_format(style='plain')
     136            p.setp(ax.get_xticklabels(), visible=False)
     137            p.setp(ax.get_yticklabels(), visible=False)
     138    # }}}
     139    # {{{ box
     140    if options.exist('box'):
     141        eval(options.getfieldvalue('box'))
     142    # }}}
     143    # {{{ xlim, ylim, zlim
     144    if options.exist('xlim'):
     145        ax.set_xlim(options.getfieldvalue('xlim'))
     146    if options.exist('ylim'):
     147        ax.set_ylim(options.getfieldvalue('ylim'))
     148    if options.exist('zlim'):
     149        ax.set_zlim(options.getfieldvalue('zlim'))
     150    # }}}
     151    # {{{ latlon TODO
     152    # }}}
     153    # {{{ Basinzoom TODO
     154    # }}}
     155    # {{{ ShowBasins TODO
     156    # }}}
     157    # {{{ clim
     158    if options.exist('clim'):
     159        lims=options.getfieldvalue('clim')
     160        assert len(lims)==2, 'error, clim should be passed as a list of length 2'
     161    elif options.exist('caxis'):
     162        lims=options.getfieldvalue('caxis')
     163        assert len(lims)==2, 'error, caxis should be passed as a list of length 2'
     164        options.addfielddefault('clim',lims)
     165    else:
     166        if len(data)>0:
     167            lims=[data.min(),data.max()]
     168        else:
     169            lims=[0,1]
     170    # }}}
     171    # {{{ shading TODO
     172    #if options.exist('shading'):
     173    # }}}
     174    # {{{ grid
     175    if options.exist('grid'):
     176        if 'on' in options.getfieldvalue('grid','on'):
     177            ax.grid()
     178    # }}}
     179    # {{{ colormap
     180    if options.exist('colornorm'):
     181        norm=options.getfieldvalue('colornorm')
     182    if options.exist('colormap'):
     183        cmap=getcolormap(options)
     184    cbar_extend=0
     185    if options.exist('cmap_set_over'):
     186        cbar_extend+=1
     187    if options.exist('cmap_set_under'):
     188        cbar_extend+=2
     189    # }}}
     190    # {{{ contours
     191    if options.exist('contourlevels'):
     192        plot_contour(md,data,options,ax)
     193    # }}}
     194    # {{{ wrapping TODO
     195    # }}}
     196    # {{{ colorbar
     197    if options.getfieldvalue('colorbar',1)==1:
     198        if cbar_extend==0:
     199            extend='neither'
     200        elif cbar_extend==1:
     201            extend='max'
     202        elif cbar_extend==2:
     203            extend='min'
     204        elif cbar_extend==3:
     205            extend='both'
     206
     207        cb = mpl.colorbar.ColorbarBase(ax.cax,cmap=cmap, norm=norm, extend=extend)
     208        if options.exist('alpha'):
     209            cb.set_alpha(options.getfieldvalue('alpha'))
     210        if options.exist('colorbarnumticks'):
     211            cb.locator=MaxNLocator(nbins=options.getfieldvalue('colorbarnumticks',5))
     212        else:
     213            cb.locator=MaxNLocator(nbins=5) # default 5 ticks
     214        if options.exist('colorbartickspacing'):
     215            locs=np.arange(lims[0],lims[1]+1,options.getfieldvalue('colorbartickspacing'))
     216            cb.set_ticks(locs)
     217        if options.exist('colorbarlines'):
     218            locs=np.arange(lims[0],lims[1]+1,options.getfieldvalue('colorbarlines'))
     219            cb.add_lines(locs,['k' for i in range(len(locs))],np.ones_like(locs))
     220        if options.exist('colorbarlineatvalue'):
     221            locs=options.getfieldvalue('colorbarlineatvalue')
     222            colors=options.getfieldvalue('colorbarlineatvaluecolor',['k' for i in range (len(locs))])
     223            widths=options.getfieldvalue('colorbarlineatvaluewidth',np.ones_like(locs))
     224            cb.add_lines(locs,colors,widths)
     225        if options.exist('colorbartitle'):
     226            if options.exist('colorbartitlepad'):
     227                cb.set_label(options.getfieldvalue('colorbartitle'),
     228                                         labelpad=options.getfieldvalue('colorbartitlepad'),fontsize=fontsize)
     229            else:
     230                cb.set_label(options.getfieldvalue('colorbartitle'),fontsize=fontsize)
     231        cb.ax.tick_params(labelsize=fontsize)
     232        cb.solids.set_rasterized(True)
     233        cb.update_ticks()
     234        cb.draw_all()
     235        if options.exist('colorbarfontsize'):
     236            colorbarfontsize=options.getfieldvalue('colorbarfontsize')
     237            cb.ax.tick_params(labelsize=colorbarfontsize)
     238            # cb.set_ticks([0,-10])
     239            # cb.set_ticklabels([-10,0,10])
     240        if options.exist('colorbarticks'):
     241            colorbarticks=options.getfieldvalue('colorbarticks')
     242            cb.set_ticks(colorbarticks)
     243        plt.sca(ax) # return to original axes control
     244    # }}}
     245    # {{{ expdisp
     246    if options.exist('expdisp'):
     247         expdisp(ax,options)
     248    # }}}
     249    # {{{ area TODO
     250    # }}}
     251    # {{{ text
     252    if options.exist('text'):
     253        text=options.getfieldvalue('text')
     254        textx=options.getfieldvalue('textx')
     255        texty=options.getfieldvalue('texty')
     256        textcolor=options.getfieldvalue('textcolor')
     257        textweight=options.getfieldvalue('textweight')
     258        textrotation=options.getfieldvalue('textrotation')
     259        textfontsize=options.getfieldvalue('textfontsize')
     260        for label,x,y,size,color,weight,rotation in zip(text,textx,texty,textfontsize,textcolor,textweight,textrotation):
     261            ax.text(x,y,label,transform=ax.transAxes,fontsize=size,color=color,weight=weight,rotation=rotation)
     262    # }}}
     263    # {{{ north arrow TODO
     264    # }}}
     265    # {{{ scale ruler TODO
     266    # }}}
     267    # {{{ streamlines TOFIX
     268    if options.exist('streamlines'):
     269        plot_streamlines(md,options,ax)
     270    # }}}
     271    # {{{ axis positions TODO
     272    # }}}
     273    # {{{ figure position TODO
     274    # }}}
     275    # {{{ axes position TODO
     276    # }}}
     277    # {{{ showregion TODO
     278    # }}}
     279    # {{{ flat edges of a partition TODO
     280    # }}}
     281    # {{{ scatter TODO
     282    # }}}
     283    # {{{ backgroundcolor TODO
     284    # }}}
     285    # {{{ figurebackgroundcolor TODO
     286    # }}}
     287    # {{{ lighting TODO
     288    # }}}
     289    # {{{ point cloud TODO
     290    # }}}
     291    # {{{ inset TODO
     292    # }}}

  • issm/trunk-jpl/src/m/plot/colormaps/cmaptools.py

    r23716 r23920  
    11import numpy as  np
     2from demmap import *
    23
    34try:
    45        import matplotlib as mpl
     6        import matplotlib.pyplot as plt
     7        from matplotlib.colors import ListedColormap, rgb_to_hsv, hsv_to_rgb
    58except ImportError:
    69        print('cannot import matplotlib, no plotting capabilities enabled')
     10
     11def getcolormap(options):
     12    '''
     13    get colormap from options and apply
     14   
     15    default: viridis
     16    supported:
     17        matplotlib defaults (see: pyplot.colormaps())
     18        Rignot
     19        demmap(50,-300,1200)
     20        demmap(50,-300,1200,'ibcao')
     21   
     22    Usage:
     23        cmap = getcolormap(options)
     24    '''
     25   
     26    map_name = options.getfieldvalue('colormap')
     27    cmap = 'viridis'
     28   
     29    # already a valid colormap, the name of a valid colormap, or empty (use default)
     30    if type(map_name) == mpl.colors.ListedColormap:
     31        return map_name
     32    elif map_name in plt.colormaps():
     33        return map_name
     34    elif map_name == '':
     35        return cmap
     36   
     37    # if we don't have a matching colormap, build one
     38    if map_name == 'Rignot':
     39        alpha=options.getfieldvalue('alpha',1)
     40        cmap = np.array((np.linspace(0,1,128,False),np.ones(128,),np.ones(128,))).T
     41        cmap[:,1] =np.maximum(np.minimum((0.1+cmap[:,0]**(1/alpha)),1),0)
     42        cmap = hsv_to_rgb(cmap)
     43        # construct a colormap object from an array of shape (n,3/4)
     44        cmap = ListedColormap(cmap)
     45       
     46    #elif map_name == 'Ala':
     47       
     48    else:
     49        # map is a library or executable function that constructs a colormap,
     50        #   function must be imported above
     51        try:
     52            cmap = ListedColormap(eval(map_name))
     53        except:
     54            raise RuntimeError("getcolormap: Error: provided colormap must be supported map or map-constructing function with syntax: 'jet' or 'function(args)'")
     55
     56    return cmap
     57
    758
    859def truncate_colormap(cmap, minval=0.0, maxval=1.0, n=100):

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

    r23716 r23920  
    1 from cmaptools import truncate_colormap
     1from cmaptools import getcolormap,truncate_colormap
    22from plot_quiver import plot_quiver
    33from scipy.interpolate import griddata
    44import numpy as  np
    55try:
    6         import matplotlib as mpl
    7         import matplotlib.pyplot as plt
    8         from mpl_toolkits.axes_grid1 import inset_locator
    9         from mpl_toolkits.mplot3d import Axes3D
    10         from mpl_toolkits.mplot3d.art3d import Poly3DCollection
     6    import matplotlib as mpl
     7    import matplotlib.pyplot as plt
     8    from mpl_toolkits.axes_grid1 import inset_locator
     9    from mpl_toolkits.mplot3d import Axes3D
     10    from mpl_toolkits.mplot3d.art3d import Poly3DCollection
    1111except ImportError:
    12         print("could not import pylab, matplotlib has not been installed, no plotting capabilities enabled")
     12    print("could not import pylab, matplotlib has not been installed, no plotting capabilities enabled")
    1313
    1414def plot_unit(x,y,z,elements,data,is2d,isplanet,datatype,options,fig,axgrid,gridindex):
    15         """
    16         PLOT_UNIT - unit plot, display data
    17 
    18         Usage:
    19         plot_unit(x,y,z,elements,data,is2d,isplanet,datatype,options)
    20 
    21         See also: PLOTMODEL, PLOT_MANAGER
    22         """
    23         #if we are plotting 3d replace the current axis
    24         if not is2d:
    25                 axgrid[gridindex].axis('off')
    26                 ax=inset_locator.inset_axes(axgrid[gridindex],width='100%',height='100%',loc=3,borderpad=0,axes_class=Axes3D)
    27                 ax.set_axis_bgcolor((0.7,0.7,0.7))
    28         else:
    29                 ax=axgrid[gridindex]
    30 
    31         #edgecolor
    32         edgecolor=options.getfieldvalue('edgecolor','None')
    33 
    34         # colormap
    35         # {{{ give number of colorlevels and transparency
    36         colorlevels=options.getfieldvalue('colorlevels',128)
    37         alpha=options.getfieldvalue('alpha',1)
    38         if alpha<1:
    39                 antialiased=True
    40         else:
    41                 antialiased=False
    42         # }}}
    43         # {{{ define wich colormap to use
    44         try:
    45                 defaultmap=plt.cm.get_cmap('viridis',colorlevels)
    46         except AttributeError:
    47                 print("Viridis can't be found (probably too old Matplotlib) reverting to gnuplot colormap")
    48                 defaultmap=truncate_colormap('gnuplot2',0.1,0.9,colorlevels)
    49         cmap=options.getfieldvalue('colormap',defaultmap)
    50         if options.exist('cmap_set_over'):
    51                 over=options.getfieldvalue('cmap_set_over','0.5')
    52                 cmap.set_over(over)
    53         if options.exist('cmap_set_under'):
    54                 under=options.getfieldvalue('cmap_set_under','0.5')
    55                 cmap.set_under(under)
    56         options.addfield('colormap',cmap)
    57         # }}}
    58         # {{{ if plotting only one of several layers reduce dataset, same for surface
    59         if options.getfieldvalue('layer',0)>=1:
    60                 plotlayer=options.getfieldvalue('layer',0)
    61                 if datatype==1:
    62                         slicesize=np.shape(elements)[0]
    63                 elif datatype in [2,3]:
    64                         slicesize=len(x)
    65                 data=data[(plotlayer-1)*slicesize:plotlayer*slicesize]
    66         # }}}
    67         # {{{ Get the colormap limits
    68         if options.exist('clim'):
    69                 lims=options.getfieldvalue('clim',[np.amin(data),np.amax(data)])
    70         elif options.exist('caxis'):
    71                 lims=options.getfieldvalue('caxis',[np.amin(data),np.amax(data)])
    72         else:
    73                 if np.amin(data)==np.amax(data):
    74                         lims=[np.amin(data)-0.5,np.amax(data)+0.5]
    75                 else:
    76                         lims=[np.amin(data),np.amax(data)]
    77         # }}}
    78         # {{{ Set the spread of the colormap (default is normal
    79         if options.exist('log'):
    80                 norm = mpl.colors.LogNorm(vmin=lims[0], vmax=lims[1])
    81         else:
    82                 norm = mpl.colors.Normalize(vmin=lims[0], vmax=lims[1])
    83         if options.exist('log'):
    84                 norm = mpl.colors.LogNorm(vmin=lims[0], vmax=lims[1])
    85         else:
    86                 norm = mpl.colors.Normalize(vmin=lims[0], vmax=lims[1])
    87         options.addfield('colornorm',norm)
    88         # }}}
    89 
    90         # Plot depending on the datatype
    91         # {{{ data are on elements
    92 
    93         if datatype==1:
    94                 if is2d:
    95                         if options.exist('mask'):
    96                                 triangles=mpl.tri.Triangulation(x,y,elements,data.mask)
    97                         else:
    98                                 triangles=mpl.tri.Triangulation(x,y,elements)
    99                         tri=ax.tripcolor(triangles,data,colorlevels,cmap=cmap,norm=norm,alpha=alpha,edgecolors=edgecolor)
    100                 else:
    101                         #first deal with colormap
    102                         loccmap = plt.cm.ScalarMappable(cmap=cmap)
    103                         loccmap.set_array([min(data),max(data)])
    104                         loccmap.set_clim(vmin=min(data),vmax=max(data))
    105 
    106                         #dealing with prism sides
    107                         recface=np.vstack((elements[:,0],elements[:,1],elements[:,4],elements[:,3])).T
    108                         eltind=np.arange(0,np.shape(elements)[0])
    109                         recface=np.vstack((recface,np.vstack((elements[:,1],elements[:,2],elements[:,5],elements[:,4])).T))
    110                         eltind=np.hstack((eltind,np.arange(0,np.shape(elements)[0])))
    111                         recface=np.vstack((recface,np.vstack((elements[:,2],elements[:,0],elements[:,3],elements[:,5])).T))
    112                         eltind=np.hstack((eltind,np.arange(0,np.shape(elements)[0])))
    113                         tmp = np.ascontiguousarray(np.sort(recface)).view(np.dtype((np.void, recface.dtype.itemsize * recface.shape[1])))
    114                         _, idx, recur = np.unique(tmp, return_index=True, return_counts=True)
    115                         recel= recface[idx[np.where(recur==1)]]
    116                         recindex=eltind[idx[np.where(recur==1)]]
    117                         for i,rectangle in enumerate(recel):
    118                                 rec=list(zip(x[rectangle],y[rectangle],z[rectangle]))
    119                                 pl3=Poly3DCollection([rec])
    120                                 color=loccmap.to_rgba(data[recindex[i]])
    121                                 pl3.set_edgecolor(color)
    122                                 pl3.set_color(color)
    123                                 ax.add_collection3d(pl3)
    124 
    125                         #dealing with prism bases
    126                         triface=np.vstack((elements[:,0:3],elements[:,3:6]))
    127                         eltind=np.arange(0,np.shape(elements)[0])
    128                         eltind=np.hstack((eltind,np.arange(0,np.shape(elements)[0])))
    129                         tmp = np.ascontiguousarray(triface).view(np.dtype((np.void, triface.dtype.itemsize * triface.shape[1])))
    130                         _, idx,recur = np.unique(tmp, return_index=True,return_counts=True)
    131                         #we keep only top and bottom elements
    132                         triel= triface[idx[np.where(recur==1)]]
    133                         triindex=eltind[idx[np.where(recur==1)]]
    134                         for i,triangle in enumerate(triel):
    135                                 tri=list(zip(x[triangle],y[triangle],z[triangle]))
    136                                 pl3=Poly3DCollection([tri])
    137                                 color=loccmap.to_rgba(data[triindex[i]])
    138                                 pl3.set_edgecolor(color)
    139                                 pl3.set_color(color)
    140                                 ax.add_collection3d(pl3)
    141        
    142                         ax.set_xlim([min(x),max(x)])
    143                         ax.set_ylim([min(y),max(y)])
    144                         ax.set_zlim([min(z),max(z)])
    145 
    146                         #raise ValueError('plot_unit error: 3D element plot not supported yet')
    147                 return
    148 
    149         # }}}
    150         # {{{ data are on nodes
    151 
    152         elif datatype==2:
    153                 if is2d:
    154                         if np.ma.is_masked(data):
    155                                 EltMask=np.asarray([np.any(np.in1d(index,np.where(data.mask))) for index in elements])
    156                                 triangles=mpl.tri.Triangulation(x,y,elements,EltMask)
    157                         else:
    158                                 triangles=mpl.tri.Triangulation(x,y,elements)
    159                                 #tri=ax.tricontourf(triangles,data,colorlevels,cmap=cmap,norm=norm,alpha=alpha)
    160                         if options.exist('log'):
    161                                 if alpha<1:#help with antialiasing
    162                                         tri=ax.tricontour(triangles,data,colorlevels,cmap=cmap,norm=norm,alpha=0.1,antialiased=antialiased)
    163                                 tri=ax.tricontourf(triangles,data,colorlevels,cmap=cmap,norm=norm,alpha=alpha,antialiased=antialiased)
    164                         else:
    165                                 if alpha<1:#help with antialiasing
    166                                         tri=ax.tricontour(triangles,data,colorlevels,cmap=cmap,norm=norm,alpha=0.1,antialiased=antialiased)
    167                                 tri=ax.tricontourf(triangles,data,colorlevels,cmap=cmap,norm=norm,alpha=alpha,extend='both',antialiased=antialiased)
    168                         if edgecolor != 'None':
    169                                 ax.triplot(x,y,elements,color=edgecolor)
    170                 else:
    171                         #first deal with the colormap
    172                         loccmap = plt.cm.ScalarMappable(cmap=cmap)
    173                         loccmap.set_array([min(data),max(data)])
    174                         loccmap.set_clim(vmin=min(data),vmax=max(data))
    175 
    176                         #deal with prism sides
    177                         recface=np.vstack((elements[:,0],elements[:,1],elements[:,4],elements[:,3])).T
    178                         recface=np.vstack((recface,np.vstack((elements[:,1],elements[:,2],elements[:,5],elements[:,4])).T))
    179                         recface=np.vstack((recface,np.vstack((elements[:,2],elements[:,0],elements[:,3],elements[:,5])).T))
    180                         tmp = np.ascontiguousarray(np.sort(recface)).view(np.dtype((np.void, recface.dtype.itemsize * recface.shape[1])))
    181                         _, idx, recur = np.unique(tmp, return_index=True, return_counts=True)
    182                         recel= recface[idx[np.where(recur==1)]]
    183                         for rectangle in recel:
    184                                 rec=list(zip(x[rectangle],y[rectangle],z[rectangle]))
    185                                 pl3=Poly3DCollection([rec])
    186                                 color=loccmap.to_rgba(np.mean(data[rectangle]))
    187                                 pl3.set_edgecolor(color)
    188                                 pl3.set_color(color)
    189                                 ax.add_collection3d(pl3)
    190 
    191                         #deal with prism faces
    192                         triface=np.vstack((elements[:,0:3],elements[:,3:6]))
    193                         tmp = np.ascontiguousarray(triface).view(np.dtype((np.void, triface.dtype.itemsize * triface.shape[1])))
    194                         _, idx,recur = np.unique(tmp, return_index=True,return_counts=True)
    195                         #we keep only top and bottom elements
    196                         triel= triface[idx[np.where(recur==1)]]
    197                         for triangle in triel:
    198                                 tri=list(zip(x[triangle],y[triangle],z[triangle]))
    199                                 pl3=Poly3DCollection([tri])
    200                                 color=loccmap.to_rgba(np.mean(data[triangle]))
    201                                 pl3.set_edgecolor(color)
    202                                 pl3.set_color(color)
    203                                 ax.add_collection3d(pl3)
    204 
    205                         ax.set_xlim([min(x),max(x)])
    206                         ax.set_ylim([min(y),max(y)])
    207                         ax.set_zlim([min(z),max(z)])
    208                         #raise ValueError('plot_unit error: 3D element plot not supported yet')
    209                 return
    210 
    211         # }}}
    212         # {{{ plotting quiver
    213         elif datatype==3:
    214                 if is2d:
    215                         Q=plot_quiver(x,y,data,options,ax)
    216                 else:
    217                         raise ValueError('plot_unit error: 3D node plot not supported yet')
    218                 return
    219        
    220         # }}}
    221         # {{{ plotting P1 Patch (TODO)
    222 
    223         elif datatype==4:
    224                 print('plot_unit message: P1 patch plot not implemented yet')
    225                 return
    226 
    227         # }}}
    228         # {{{ plotting P0 Patch (TODO)
    229 
    230         elif datatype==5:
    231                 print('plot_unit message: P0 patch plot not implemented yet')
    232                 return
    233 
    234         # }}}
    235         else:
    236                 raise ValueError('datatype=%d not supported' % datatype)
     15    """
     16    PLOT_UNIT - unit plot, display data
     17
     18    Usage:
     19    plot_unit(x,y,z,elements,data,is2d,isplanet,datatype,options)
     20
     21    See also: PLOTMODEL, PLOT_MANAGER
     22    """
     23    #if we are plotting 3d replace the current axis
     24    if not is2d:
     25        axgrid[gridindex].axis('off')
     26        ax=inset_locator.inset_axes(axgrid[gridindex],width='100%',height='100%',loc=3,borderpad=0,axes_class=Axes3D)
     27        ax.set_axis_bgcolor((0.7,0.7,0.7))
     28    else:
     29        ax=axgrid[gridindex]
     30
     31    #edgecolor
     32    edgecolor=options.getfieldvalue('edgecolor','None')
     33
     34    # colormap
     35    # {{{ give number of colorlevels and transparency
     36    colorlevels=options.getfieldvalue('colorlevels',128)
     37    alpha=options.getfieldvalue('alpha',1)
     38    if alpha<1:
     39        antialiased=True
     40    else:
     41        antialiased=False
     42    # }}}
     43    # {{{ define wich colormap to use
     44    try:
     45        defaultmap=plt.cm.get_cmap('viridis',colorlevels)
     46    except AttributeError:
     47        print("Viridis can't be found (probably too old Matplotlib) reverting to gnuplot colormap")
     48        defaultmap=truncate_colormap('gnuplot2',0.1,0.9,colorlevels)
     49    if not options.exist('colormap'):
     50        cmap=defaultmap
     51    else:
     52        cmap=getcolormap(options)
     53    if options.exist('cmap_set_over'):
     54        over=options.getfieldvalue('cmap_set_over','0.5')
     55        cmap.set_over(over)
     56    if options.exist('cmap_set_under'):
     57        under=options.getfieldvalue('cmap_set_under','0.5')
     58        cmap.set_under(under)
     59    options.addfield('colormap',cmap)
     60    # }}}
     61    # {{{ if plotting only one of several layers reduce dataset, same for surface
     62    if options.getfieldvalue('layer',0)>=1:
     63        plotlayer=options.getfieldvalue('layer',0)
     64        if datatype==1:
     65            slicesize=np.shape(elements)[0]
     66        elif datatype in [2,3]:
     67            slicesize=len(x)
     68        data=data[(plotlayer-1)*slicesize:plotlayer*slicesize]
     69    # }}}
     70    # {{{ Get the colormap limits
     71    if options.exist('clim'):
     72        lims=options.getfieldvalue('clim',[np.amin(data),np.amax(data)])
     73    elif options.exist('caxis'):
     74        lims=options.getfieldvalue('caxis',[np.amin(data),np.amax(data)])
     75    else:
     76        if np.amin(data)==np.amax(data):
     77            lims=[np.amin(data)-0.5,np.amax(data)+0.5]
     78        else:
     79            lims=[np.amin(data),np.amax(data)]
     80    # }}}
     81    # {{{ Set the spread of the colormap (default is normal
     82    if options.exist('log'):
     83        norm = mpl.colors.LogNorm(vmin=lims[0], vmax=lims[1])
     84    else:
     85        norm = mpl.colors.Normalize(vmin=lims[0], vmax=lims[1])
     86    if options.exist('log'):
     87        norm = mpl.colors.LogNorm(vmin=lims[0], vmax=lims[1])
     88    else:
     89        norm = mpl.colors.Normalize(vmin=lims[0], vmax=lims[1])
     90    options.addfield('colornorm',norm)
     91    # }}}
     92
     93    # Plot depending on the datatype
     94    # {{{ data are on elements
     95
     96    if datatype==1:
     97        if is2d:
     98            if options.exist('mask'):
     99                triangles=mpl.tri.Triangulation(x,y,elements,data.mask)
     100            else:
     101                triangles=mpl.tri.Triangulation(x,y,elements)
     102            tri=ax.tripcolor(triangles,data,colorlevels,cmap=cmap,norm=norm,alpha=alpha,edgecolors=edgecolor)
     103        else:
     104            #first deal with colormap
     105            loccmap = plt.cm.ScalarMappable(cmap=cmap)
     106            loccmap.set_array([min(data),max(data)])
     107            loccmap.set_clim(vmin=min(data),vmax=max(data))
     108
     109            #dealing with prism sides
     110            recface=np.vstack((elements[:,0],elements[:,1],elements[:,4],elements[:,3])).T
     111            eltind=np.arange(0,np.shape(elements)[0])
     112            recface=np.vstack((recface,np.vstack((elements[:,1],elements[:,2],elements[:,5],elements[:,4])).T))
     113            eltind=np.hstack((eltind,np.arange(0,np.shape(elements)[0])))
     114            recface=np.vstack((recface,np.vstack((elements[:,2],elements[:,0],elements[:,3],elements[:,5])).T))
     115            eltind=np.hstack((eltind,np.arange(0,np.shape(elements)[0])))
     116            tmp = np.ascontiguousarray(np.sort(recface)).view(np.dtype((np.void, recface.dtype.itemsize * recface.shape[1])))
     117            _, idx, recur = np.unique(tmp, return_index=True, return_counts=True)
     118            recel= recface[idx[np.where(recur==1)]]
     119            recindex=eltind[idx[np.where(recur==1)]]
     120            for i,rectangle in enumerate(recel):
     121                rec=list(zip(x[rectangle],y[rectangle],z[rectangle]))
     122                pl3=Poly3DCollection([rec])
     123                color=loccmap.to_rgba(data[recindex[i]])
     124                pl3.set_edgecolor(color)
     125                pl3.set_color(color)
     126                ax.add_collection3d(pl3)
     127
     128            #dealing with prism bases
     129            triface=np.vstack((elements[:,0:3],elements[:,3:6]))
     130            eltind=np.arange(0,np.shape(elements)[0])
     131            eltind=np.hstack((eltind,np.arange(0,np.shape(elements)[0])))
     132            tmp = np.ascontiguousarray(triface).view(np.dtype((np.void, triface.dtype.itemsize * triface.shape[1])))
     133            _, idx,recur = np.unique(tmp, return_index=True,return_counts=True)
     134            #we keep only top and bottom elements
     135            triel= triface[idx[np.where(recur==1)]]
     136            triindex=eltind[idx[np.where(recur==1)]]
     137            for i,triangle in enumerate(triel):
     138                tri=list(zip(x[triangle],y[triangle],z[triangle]))
     139                pl3=Poly3DCollection([tri])
     140                color=loccmap.to_rgba(data[triindex[i]])
     141                pl3.set_edgecolor(color)
     142                pl3.set_color(color)
     143                ax.add_collection3d(pl3)
     144   
     145            ax.set_xlim([min(x),max(x)])
     146            ax.set_ylim([min(y),max(y)])
     147            ax.set_zlim([min(z),max(z)])
     148
     149            #raise ValueError('plot_unit error: 3D element plot not supported yet')
     150        return
     151
     152    # }}}
     153    # {{{ data are on nodes
     154
     155    elif datatype==2:
     156        if is2d:
     157            if np.ma.is_masked(data):
     158                EltMask=np.asarray([np.any(np.in1d(index,np.where(data.mask))) for index in elements])
     159                triangles=mpl.tri.Triangulation(x,y,elements,EltMask)
     160            else:
     161                triangles=mpl.tri.Triangulation(x,y,elements)
     162                #tri=ax.tricontourf(triangles,data,colorlevels,cmap=cmap,norm=norm,alpha=alpha)
     163            if options.exist('log'):
     164                if alpha<1:#help with antialiasing
     165                    tri=ax.tricontour(triangles,data,colorlevels,cmap=cmap,norm=norm,alpha=0.1,antialiased=antialiased)
     166                tri=ax.tricontourf(triangles,data,colorlevels,cmap=cmap,norm=norm,alpha=alpha,antialiased=antialiased)
     167            else:
     168                if alpha<1:#help with antialiasing
     169                    tri=ax.tricontour(triangles,data,colorlevels,cmap=cmap,norm=norm,alpha=0.1,antialiased=antialiased)
     170                tri=ax.tricontourf(triangles,data,colorlevels,cmap=cmap,norm=norm,alpha=alpha,extend='both',antialiased=antialiased)
     171            if edgecolor != 'None':
     172                ax.triplot(x,y,elements,color=edgecolor)
     173        else:
     174            #first deal with the colormap
     175            loccmap = plt.cm.ScalarMappable(cmap=cmap)
     176            loccmap.set_array([min(data),max(data)])
     177            loccmap.set_clim(vmin=min(data),vmax=max(data))
     178
     179            #deal with prism sides
     180            recface=np.vstack((elements[:,0],elements[:,1],elements[:,4],elements[:,3])).T
     181            recface=np.vstack((recface,np.vstack((elements[:,1],elements[:,2],elements[:,5],elements[:,4])).T))
     182            recface=np.vstack((recface,np.vstack((elements[:,2],elements[:,0],elements[:,3],elements[:,5])).T))
     183            tmp = np.ascontiguousarray(np.sort(recface)).view(np.dtype((np.void, recface.dtype.itemsize * recface.shape[1])))
     184            _, idx, recur = np.unique(tmp, return_index=True, return_counts=True)
     185            recel= recface[idx[np.where(recur==1)]]
     186            for rectangle in recel:
     187                rec=list(zip(x[rectangle],y[rectangle],z[rectangle]))
     188                pl3=Poly3DCollection([rec])
     189                color=loccmap.to_rgba(np.mean(data[rectangle]))
     190                pl3.set_edgecolor(color)
     191                pl3.set_color(color)
     192                ax.add_collection3d(pl3)
     193
     194            #deal with prism faces
     195            triface=np.vstack((elements[:,0:3],elements[:,3:6]))
     196            tmp = np.ascontiguousarray(triface).view(np.dtype((np.void, triface.dtype.itemsize * triface.shape[1])))
     197            _, idx,recur = np.unique(tmp, return_index=True,return_counts=True)
     198            #we keep only top and bottom elements
     199            triel= triface[idx[np.where(recur==1)]]
     200            for triangle in triel:
     201                tri=list(zip(x[triangle],y[triangle],z[triangle]))
     202                pl3=Poly3DCollection([tri])
     203                color=loccmap.to_rgba(np.mean(data[triangle]))
     204                pl3.set_edgecolor(color)
     205                pl3.set_color(color)
     206                ax.add_collection3d(pl3)
     207
     208            ax.set_xlim([min(x),max(x)])
     209            ax.set_ylim([min(y),max(y)])
     210            ax.set_zlim([min(z),max(z)])
     211            #raise ValueError('plot_unit error: 3D element plot not supported yet')
     212        return
     213
     214    # }}}
     215    # {{{ plotting quiver
     216    elif datatype==3:
     217        if is2d:
     218            Q=plot_quiver(x,y,data,options,ax)
     219        else:
     220            raise ValueError('plot_unit error: 3D node plot not supported yet')
     221        return
     222   
     223    # }}}
     224    # {{{ plotting P1 Patch (TODO)
     225
     226    elif datatype==4:
     227        print('plot_unit message: P1 patch plot not implemented yet')
     228        return
     229
     230    # }}}
     231    # {{{ plotting P0 Patch (TODO)
     232
     233    elif datatype==5:
     234        print('plot_unit message: P0 patch plot not implemented yet')
     235        return
     236
     237    # }}}
     238    else:
     239        raise ValueError('datatype=%d not supported' % datatype)
Note: See TracChangeset for help on using the changeset viewer.