defangle.py

   1 import matplotlib as mpl
   2 import matplotlib.mlab as mlab
   3 import matplotlib.pyplot as plt
   4 import matplotlib.axes as axe
   5 import matplotlib.collections as collections
   6 import numpy as np
   7 import random
   8 from pylab import *
   9 from matplotlib.pyplot import *
  10
  11 ANGLE = []
  12 VELO = []
  13 ANGLE.append("../rk_defangle.dat")
  14 ANGLE.append("../data/eu_defangle_t225.dat")
  15 #ANGLE.append("../data/eu_defangle_t675.dat")
  16 VELO.append("../rk_bif_avg_velo.dat")
  17 VELO.append("../data/eu_defangle_velo_t255.dat")
  18 #VELO.append("../data/eu_defangle_velo_t675.dat")
  19
  20
  21 FOUT = "../defangle.png"
  22
  23 xlabel = r'$h$'
  24 ylabel = r'$\delta y$'
  25 res = (1024/80,768/80)
  26 pad1 = 1
  27 pad2 = 0.1
  28 tsize = 16
  29 lsize = 24
  30 max = 0.875
  31 acangle = 58
  32
  33 ######################################################################
  34
  35 angle = []
  36 velo = []
  37
  38 for i in range(len(ANGLE)):
  39         angle.append(mlab.csv2rec(ANGLE[i], delimiter='\t', comments='#'))
  40 for i in range(len(VELO)):
  41         velo.append(mlab.csv2rec(VELO[i], delimiter='\t', comments='#'))
  42
  43 r = []
  44 psi = []
  45 vx = []
  46
  47 for i in range(len(angle)):
  48         r.append(angle[i].r)
  49         psi.append(angle[i].x)
  50 for i in range(len(velo)):
  51         vx.append(velo[i].x)
  52
  53
  54
  55 def mean_angle(param,angle):
  56         unique_param = unique(param)
  57         list_angles = [ [] for DUMMYVAR in range(len(unique_param)) ]
  58         ret_avg = []
  59         for i in range(len(unique_param)):
  60                 for j in range(len(angle)):
  61                         #print j, len(param), i, len(unique_param)
  62                         if (param[j] == unique_param[i]):
  63                                 list_angles[i].append(angle[j])
  64                                 #print len(list_angles[i])
  65                 ret_avg.append(abs(float(sum(list_angles[i])) /
  66                         len(list_angles[i])))
  67                 #print len(ret_avg)
  68         return ret_avg
  69
  70 marker_list = [
  71                 '.',
  72                 ',',
  73                 'o',
  74                 'v',
  75                 ]
  76
  77 ls_list = [
  78                 '-',
  79                 '--',
  80                 '-.',
  81                 ':',
  82                 '.',
  83                 ',',
  84                 'o',
  85                 'v',
  86                 ]
  87
  88 print marker_list[0]
  89
  90 left, width = 0.1, 0.8
  91 c = max/2
  92 rect1 = [left, c+0.1, width, c]
  93 rect2 = [left, 0.1, width, c]
  94
  95 rad_ac = acangle*(np.pi/180)
  96
  97 fig = plt.figure(figsize=res)
  98
  99 ax1 = fig.add_axes(rect1)
 100 ax2 = fig.add_axes(rect2)
 101
 102 avg_angles = []
 103 avg_velo = []
 104 unique_r = []
 105
 106 for i in range(len(angle)):
 107         print "avg_angles[",i,"]"
 108         avg_angles.append(mean_angle(r[i],psi[i]))
 109         print "avg_velo[",i,"]"
 110
 111         avg_velo.append(mean_angle(r[i],vx[i]))
 112         unique_r.append(unique(r[i]))
 113
 114 for i in range(len(angle)):
 115         ax1.plot(unique_r[i], avg_angles[i], marker=marker_list[i],
 116                         ms=0, ls=ls_list[i], color='black')
 117 for i in range(len(velo)):
 118         ax2.plot(unique_r[i], avg_velo[i], marker_list[i], ms=0,
 119                         ls=ls_list[i], color='black')
 120
 121 xticks = np.arange(0,360+1,45)
 122 yticks1 = np.arange(0,360+1,45)
 123 yticks2 = np.arange(0,1.1,0.5)
 124
 125 ax1.set_xlim(r[0].min()-pad1,r[0].max()+pad1)
 126 ax1.set_ylim(yticks1.min(),yticks1.max())
 127 ax2.set_xlim(r[0].min()-pad2,r[0].max()+pad2)
 128 ax2.set_ylim(0,vx[i].max()+pad2)
 129 ax1.set_xticks([])
 130 ax2.set_xticks(xticks)
 131 ax1.set_yticks(yticks1)
 132 ax2.set_yticks(yticks2)
 133 ax1.axvline(acangle, color='black', ls=":")
 134 #ax1.axvline(acangle-45, color='black', ls="-.")
 135 #ax1.axvline(acangle+45, color='black', ls="-.")
 136 ax1.axvline(acangle+180, color='black', ls=":")
 137 #ax1.axvline(acangle+180-45, color='black', ls="-.")
 138 #ax1.axvline(acangle+180+45, color='black', ls="-.")
 139 ax2.axvline(acangle, color='black', ls=":")
 140 ax2.axvline(acangle+180, color='black', ls=":")
 141 #ax2.axvline(acangle+180-45, color='black', ls="-.")
 142 #ax2.axvline(acangle+180+45, color='black', ls="-.")
 143 #ax2.axvline(acangle-45, color='black', ls="-.")
 144 #ax2.axvline(acangle+45, color='black', ls="-.")
 145
 146 ax1.set_ylabel(r'$\psi$', size=lsize)
 147 ax2.set_ylabel(r'$v_x$', size=lsize)
 148 ax2.set_xlabel(r'$\phi$', size=lsize)
 149 ax2.set_xticklabels(ax2.get_xticks(), size=tsize)
 150 ax1.set_yticklabels(ax1.get_yticks(), size=tsize)
 151 ax2.set_yticklabels(ax2.get_yticks(), size=tsize)
 152
 153 plt.savefig(FOUT)
 154
 155 #plt.show()