|
version 1.1, 2018/12/09 13:52:19
|
version 1.2, 2019/02/22 22:26:26
|
|
|
|
| try: | try: |
| freqscale | freqscale |
| except: | except: |
| freqscale=1.0/3000 |
freqscale=1.0/3 |
| | |
| try: | try: |
| dpi | dpi |
|
|
|
| nframes = 64 | nframes = 64 |
| | |
| try: | try: |
| animateflag |
colorshift |
| except: | except: |
| animateflag=1 |
colorshift=0 |
| |
|
| |
try: |
| |
animate |
| |
except: |
| |
animate=1 |
| |
|
| |
try: |
| |
show |
| |
except: |
| |
show=1 |
| | |
| try: | try: |
| rsurf | rsurf |
|
|
|
| pixels=1000 | pixels=1000 |
| | |
| try: | try: |
| bangle |
|
| except: |
|
| bangle=-30 |
|
| |
|
| try: |
|
| figsize | figsize |
| except: | except: |
| figsize=5 | figsize=5 |
| | |
| # not needed if using table of surface values |
sosh.nframes=nframes |
| |
sosh.dpi=dpi |
| |
sosh.icolshift=colorshift |
| |
|
| sosh.loadmodel() | sosh.loadmodel() |
| | |
| #(phi,theta)=sosh.image2sphere(xpixels=pixels,ypixels=pixels,bangle=bangle) | #(phi,theta)=sosh.image2sphere(xpixels=pixels,ypixels=pixels,bangle=bangle) |
| (r, theta) = sosh.image2rtheta(xpixels=pixels,ypixels=pixels) | (r, theta) = sosh.image2rtheta(xpixels=pixels,ypixels=pixels) |
| x=np.cos(theta) | x=np.cos(theta) |
| |
phi=0.0*theta |
| (nx,ny)=x.shape | (nx,ny)=x.shape |
| if (rsurf < 1.0): | if (rsurf < 1.0): |
| newind = (r <= rsurf) | newind = (r <= rsurf) |
| r.mask = np.logical_not(newindex) | r.mask = np.logical_not(newindex) |
| | |
| #modeparms=np.loadtxt('modes.average') |
#modeparms=np.loadtxt('mdi.average.modes') |
| # use model values instead of measured values | # use model values instead of measured values |
| modeparms=np.loadtxt('modes.model.txt') |
#modeparms=np.loadtxt('model.surface.modes') |
| lmod=modeparms[:,0] |
#lmod=modeparms[:,0] |
| nmod=modeparms[:,1] |
#nmod=modeparms[:,1] |
| | |
| arrlist={} | arrlist={} |
| lmaxlist={} | lmaxlist={} |
|
|
|
| hclist=[] | hclist=[] |
| | |
| interval = 1.0/25 | interval = 1.0/25 |
| |
maxsave=0.0 |
| | |
| #Animation function to call | #Animation function to call |
| def sumanimate(i): |
def calcsumt(i): |
| |
|
| sumr=0.0 | sumr=0.0 |
| sumt=0.0 | sumt=0.0 |
| sump=0.0 | sump=0.0 |
| |
plotvar=sosh.plotvar |
| |
|
| |
if plotvar in ['Vr','Vmag','Vsq']: |
| for k in range(nmodes): | for k in range(nmodes): |
| sumr += freqlist[k]*rclist[k]*rlist[k]*np.exp(-1.0j*(2*np.pi*freqlist[k]*float(i) / np.float(nframes))) |
sumr += freqlist[k]*rclist[k]*rlist[k]*np.exp(-1.0j*2*np.pi*freqlist[k]*freqscale*float(i) / np.float(nframes)) |
| sumt += freqlist[k]*hclist[k]*tlist[k]*np.exp(-1.0j*(2*np.pi*freqlist[k]*float(i) / np.float(nframes))) |
if plotvar in ['Vt', 'Vh', 'Vmag','Vsq']: |
| sump += freqlist[k]*hclist[k]*plist[k]*np.exp(-1.0j*(2*np.pi*freqlist[k]*float(i) / np.float(nframes))) |
for k in range(nmodes): |
| vh2=(np.square(sumt.real) + np.square(sump.real)) |
sumt += freqlist[k]*hclist[k]*tlist[k]*np.exp(-1.0j*2*np.pi*freqlist[k]*freqscale*float(i) / np.float(nframes)) |
| v2=(np.square(sumr.real) + vh2) |
if plotvar in ['Vt', 'Vh', 'Vmag','Vsq']: |
| |
for k in range(nmodes): |
| |
sump += freqlist[k]*hclist[k]*plist[k]*np.exp(-1.0j*2*np.pi*freqlist[k]*freqscale*float(i) / np.float(nframes)) |
| |
|
| |
if plotvar in ['Vh', 'Vmag', 'Vsq']: |
| |
vh2=np.square(sumt.real) + np.square(sump.real) |
| |
if plotvar in ['Vmag', 'Vsq']: |
| |
v2=np.square(sumr.real) + vh2 |
| | |
| plotvar=sosh.plotvar |
|
| if (plotvar == 'Vr'): | if (plotvar == 'Vr'): |
| d=sumr.real | d=sumr.real |
| if (plotvar == 'Vt'): | if (plotvar == 'Vt'): |
| Line 108 def sumanimate(i): |
|
| Line 128 def sumanimate(i): |
|
| if (plotvar == 'Vsq'): | if (plotvar == 'Vsq'): |
| d=v2 | d=v2 |
| | |
| |
return d |
| |
|
| |
def sumanimate(i): |
| |
|
| |
d=calcsumt(i) |
| |
global maxsave |
| |
mn,mx=d.min(),d.max() |
| |
maxabs=np.abs([mn,mx]).max() |
| |
if (maxabs > maxsave): |
| |
maxsave=maxabs |
| |
# print("new max=%f"%maxabs, end=" ", flush=True) |
| |
|
| im.set_data(d) | im.set_data(d) |
| fig.canvas.draw() | fig.canvas.draw() |
| return |
return maxabs |
| | |
| sosh.animate=sumanimate | sosh.animate=sumanimate |
| |
sosh.calcimage=calcsumt |
| sosh.isave=0 | sosh.isave=0 |
| | |
| |
|
| print("You may enter 'q' to quit or 'c' to change parameters.") | print("You may enter 'q' to quit or 'c' to change parameters.") |
| instr=sosh.catchc("Enter number of modes: ",2) | instr=sosh.catchc("Enter number of modes: ",2) |
| if (instr == 'q'): | if (instr == 'q'): |
| Line 126 for i in range(nmodes): |
|
| Line 160 for i in range(nmodes): |
|
| (l,m,n)=sosh.querylmn() | (l,m,n)=sosh.querylmn() |
| if l == -1: | if l == -1: |
| sys.exit() | sys.exit() |
| ntest = nmod[l==lmod] |
ntest = sosh.modeln[l==sosh.modell] |
| while (n not in ntest): | while (n not in ntest): |
| print("That n was not modelled for that l. Modelled values are in the range %i to %i. Try again." \ | print("That n was not modelled for that l. Modelled values are in the range %i to %i. Try again." \ |
| % (ntest.min(),ntest.max())) | % (ntest.min(),ntest.max())) |
| (l,m,n)=sosh.querylmn() | (l,m,n)=sosh.querylmn() |
| ntest = nmod[l==lmod] |
ntest = sosh.modeln[l==sosh.modell] |
| | |
| signedm=m | signedm=m |
| m=np.abs(m) | m=np.abs(m) |
| Line 148 for i in range(nmodes): |
|
| Line 182 for i in range(nmodes): |
|
| print("Spherical harmonic retrieved.") | print("Spherical harmonic retrieved.") |
| | |
| if (m > 0): | if (m > 0): |
| phi=(np.pi/4)/m |
phi[:,int(nx/2):nx]=(np.pi/4)/m |
| |
phi[:,0:int(nx/2)]=(np.pi/4)/m + np.pi |
| else: | else: |
| phi=0.0 |
phi=0.0*theta |
| ylm=y*np.exp(1.0j*signedm*phi) | ylm=y*np.exp(1.0j*signedm*phi) |
| dylmt=-np.sin(theta)*dy*np.exp(1.0j*signedm*phi) | dylmt=-np.sin(theta)*dy*np.exp(1.0j*signedm*phi) |
| dylmp=1.0j*m*y*np.exp(1.0j*signedm*phi) |
dylmp=1.0j*signedm*y*np.exp(1.0j*signedm*phi) |
| | |
| rlist.append(ylm) | rlist.append(ylm) |
| tlist.append(dylmt) | tlist.append(dylmt) |
| Line 163 for i in range(nmodes): |
|
| Line 198 for i in range(nmodes): |
|
| fr=np.interp(r,sosh.rmesh,xir*np.sqrt(sosh.rho)) | fr=np.interp(r,sosh.rmesh,xir*np.sqrt(sosh.rho)) |
| fh=np.interp(r,sosh.rmesh,xih*np.sqrt(sosh.rho)) | fh=np.interp(r,sosh.rmesh,xih*np.sqrt(sosh.rho)) |
| | |
| ind = ((l==lmod) & (n==nmod)) |
idx = ((sosh.modeln == n) & (sosh.modell == l)) |
| freq = float(modeparms[ind,2])*freqscale |
freq = sosh.modelnu[idx]/1000.0 |
| # rc = float(modeparms[ind,5])/1e8 |
# ind = ((l==lmod) & (n==nmod)) |
| # hc = float(modeparms[ind,6])/1e8 |
# freq = float(modeparms[ind,2])/1000.0 |
| rclist.append(fr) | rclist.append(fr) |
| hclist.append(fh) | hclist.append(fh) |
| freqlist.append(freq) | freqlist.append(freq) |
| # print("Using vertical to horizontal ratio of %f, theoretical is %f." % (rc/hc, float(modeparms[ind,4]))) |
print("Scaled frequency = %f." % (freq*freqscale)) |
| print("Scaled frequency = %f." % (freq)) |
|
| | |
| # if using measured mode parameters, this would | # if using measured mode parameters, this would |
| # allow you to compute frequency as a function of m | # allow you to compute frequency as a function of m |
| Line 204 varlist['Vsq']=v2 |
|
| Line 238 varlist['Vsq']=v2 |
|
| while True: | while True: |
| | |
| var=varlist[sosh.plotvar] | var=varlist[sosh.plotvar] |
| |
maxsave=0.0 |
| | |
| fig,im = sosh.drawfigure(var, fsize=figsize) | fig,im = sosh.drawfigure(var, fsize=figsize) |
| if (sosh.isave != 0): | if (sosh.isave != 0): |
| sosh.savefigure(animateflag=animateflag) |
sosh.savefigure(ianimate=animate, label='_interior') |
| if (animateflag != 0): |
if (show != 0): |
| |
if (animate != 0): |
| anim = animation.FuncAnimation(fig, sumanimate, frames=nframes, interval=interval) | anim = animation.FuncAnimation(fig, sumanimate, frames=nframes, interval=interval) |
| plt.show() | plt.show() |
| | |
Return to
addradial.py
CVS log
Up to [Development] / JSOC / proj / globalhs / sosh