c c------------------------------------------------------------------ c Begin file /home/sting/foster/tigcm/dynglb.f c------------------------------------------------------------------ c subroutine dynglb(utc) c c Make global (cylindrical equidistant) contours c include 'tgcmparam.h' include 'input.h' include 'tigcm.h' include 'tgcmhdr.h' include 'tiegcmfld.h' include 'tigcmfld.h' include 'selgrid.h' include 'tgcmlab.h' include 'cool.h' include 'plt.h' include 'color.h' dimension plotp(imx,jmx) character*56 toplab,lab,lab1 c c Declarations for axes labeling: c parameter (nx=7,ny=7) dimension numx(nx),numy(ny) data numx/-180,-120,-60,0,60,120,180/ data numy/-90,-60,-30,0,30,60,90/ character*8 fmtx,fmty character*16 labx,laby data fmtx/'(i4) '/, fmty/'(i3) '/ data nfmtx/4/, nfmty/3/, mnrx/3/, mnry/3/ data laby/' LATITUDE (DEG) '/ data labx/' LONGITUDE (DEG)'/ data xc/-180./,xd/180./ data yc/-87.5/,yd/87.5/ c c ispvec = special value flag for velvct c spvec = specil value for velvct data ispvec/0/, spvec/1.e36/ c c Temp: data maps/1/, iccont/1/ c c Set up conpack: c call cpseti('SET',0) pmapl = 0.10 pmapr = 0.95 pmapb = 0.455 pmapt = 0.880 wl = -180. wr = +180. wb = -90. wt = +90. rilx = 0.5 rily = 0.08 rily = -((pmapb-rily) / (pmapt-pmapb)) call cpsetr('ILX - info label x-coord',rilx) call cpsetr('ILY - info label y-coord',rily) iel = 0 c c Set up ezmap, if doing continental outlines: c if (maps.le.0) then call set(pmapl,pmapr,pmapb,pmapt,wl,wr,wb,wt,1) else proj = 'CE ' plat = 0. plon = 0. call cpseti('MAP - mapping flag',1) rot = 0. call mappos(pmapl,pmapr,pmapb,pmapt) call mapsti('EL',iel) call mapstc('OU','CO') call mapsti('DO',1) call maproj(proj,plat,plon,rot) call mapset('MA',0.,0.,0.,0.) call mapsti('G1',1) call mapsti('G2',2) call mapsti('VS',0) call mapint endif c c Field loop: c write(6,"(' ')") write(6,"('Global contours (dynamo fields): npls=',i2, + ' nhts=',i2)") npls,nhts do 100 ip=1,ntiegcmf if (iptiegcm(ip).le.0) goto 100 c c Selected pressure and heights loop: c do 110 k=1,npls+nhts kk = k-npls if (k.le.npls) then write(toplab,"('UT = ',f6.2,' ZP = ',f4.1,35x)") utc,spls(k) c c Define plot array (selected pressure): plotp(:,:) = dynpnt(:,ispls(k),:,ip) else write(toplab,"('UT = ',f6.2,' HT = ',f6.1,33x)") + utc,shts(kk) c c Define plot array (interpolate to height surface): call tgcmint(dynpnt(1,1,1,ip),pnt(1,1,1,ixz),gcmlat,jmx, + gcmlon,imx,shts(kk),1,idynlog(ip),2,plotp,imx,jmx, + 1,kmx,gcmzp,ier,1,cpspval,0) endif c c Take logs if necessary: if (idynlog(ip).gt.0.and.log10pl.gt.0) then do 125 j=1,jmx do 125 i=1,imx if (plotp(i,j).eq.cpspval.or.plotp(i,j).le.1.e-30) + then plotp(i,j) = cpspval else plotp(i,j) = alog10(plotp(i,j)) endif 125 continue write(lab,"('LOG10 ',a50)") labtiegcm(ip)(1:50) else lab = labtiegcm(ip) endif c c Contour it, with axes and axis labels: c if (iclrfill.le.0) then call bwcon(plotp,imx,imx,jmx,xc,xd,yc,yd, + cintdyn(ip),cmindyn(ip),cmaxdyn(ip), + 1,nx,numx,mnrx,labx,fmtx,nfmtx, + ny,numy,mnry,laby,fmty,nfmty) else lbar = 2 call clrcon(plotp,imx,imx,jmx,xc,xd,yc,yd, + cintdyn(ip),cmindyn(ip),cmaxdyn(ip), + 1,nx,numx,mnrx,labx,fmtx,nfmtx, + ny,numy,mnry,laby,fmty,nfmty, + lbar,pmapl,pmapr,.09,.19,ip+ntigcmf) endif c call box(0) c c Add local time axis at bottom: call sltxax(utc) c c Draw continental outlines, if desired: c (white for b/w plots, iccont for color fill) c if (maps.gt.0) then call gqplci(ierr,isplci) call gqpmci(ierr,ispmci) if (iclrfill.le.0) then call gsplci(iwhite) call gspmci(iwhite) else call gsplci(iccont) call gspmci(iccont) endif call maplot call gsplci(isplci) call gspmci(ispmci) endif c c Add labels: c call top2lab(lab,0.15,toplab,0.07,.017) call botlab(runlab(1:lnblnk2(runlab)),0.5*(pmapl+pmapr),0.25) c c Bottom history label: if (iclrfill.le.0) then write(lab1,"(7x,'TGCM History = ',4a8,2x)") + (output(ii,1),ii=1,3),output(3,2) call botlab(lab1,0.5*(pmapl+pmapr),.20) endif c c Add neutral wind vectors if doing temp: c if (ip.eq.ixpot.and.ipltuv.eq.1) then if (k.le.npls) then udat(:,:) = ui(:,ispls(k),:) vdat(:,:) = vi(:,ispls(k),:) else c c Do height interpolation of un,vn: c call tgcmint(dynpnt(1,1,1,ixui),pnt(1,1,1,ixz),gcmlat,jmx, + gcmlon,imx,shts(kk),1,igcmlog(ixui),2,udat,imx,jmx, + 1,kmx,gcmzp,ier,1,cpspval,0) call tgcmint(dynpnt(1,1,1,ixvi),pnt(1,1,1,ixz),gcmlat,jmx, + gcmlon,imx,shts(kk),1,igcmlog(ixvi),2,vdat,imx,jmx, + 1,kmx,gcmzp,ier,1,cpspval,0) endif do 95 j=1,jmx do 95 i=1,imx if (udat(i,j).eq.0..and.vdat(i,j).eq.0.) then udat(i,j) = cpspval vdat(i,j) = cpspval endif 95 continue uvmax = 0. m = 1 ispvec = 0 spvec = 1.e36 iskip = 0 call pltvec(udat,vdat,uvmax,jmx,iskip,gcmlat(1), + gcmlat(jmx),m,ispvec,spvec) endif c c Finish up: call frame iframe=iframe+1 if (k.le.npls) then write(6,"('dynglb frame ',i3,': ',a,' UT=',f4.1, + ' ZP=',f4.1)") iframe,lab(1:40),utc,spls(k) else write(6,"('dynglb frame ',i3,': ',a,' UT=',f4.1, + ' HT=',f6.1)") iframe,lab(1:40),utc,shts(kk) endif c c End pressure loop: 110 continue c c End field loop 100 continue c return end