subroutine coplotvc(ele,fname,bw,be,bs,bn, *jm,b1,b2,maxplots, *olddtm,newdtm,region,filter,elecap,ij,xvlon,xvlat,xllon,xllat, *lorvog,lorvopc,igridsec,igridsec2,fn,cptlo,cpthi,cptin) c - Updated 2019 02 15 for debugging as we tried to build the digital archive. c - Something went wrong, so this version attempts to fix that. c - Updated 2018 06 29 from the version developed for VERTCON 3 (with "filter" and two grid spacings) c - to the version used for NADCON 5.01. c---------------------------- c - Updated 2018 05 16 for 2 grid spacings c - Cannibalized from NADCON 5 for VERTCON 3 c - 2016 08 29: c Expanded the refence vector calls to be 6 decimal places, as well c as the initial -R call for S/N/W/E c and also the -B part of that call c - 2016 08 25: c .gmtdefaults4 has been changed so X_ORIGIN is equal to 0.0 c Center the plot with "-Xc" at first "psxy" call c Remove all "-JM**i+" references, and just use the actual c width (jm) that came out of the "getmapbounds" routine and c was sent here. c - 2016 07 29: Updated the reference vector call to c - have the "-N" option, so it'll plot outside the map c - 2016 07 21: Updated JM usage for new forced sizes c subroutine bwplotvc(ele,fname,bw,be,bs,bn,jm,b1,b2,maxplots, c *olddtm,newdtm,region,elecap,ij,xvlon,xvlat,xllon,xllat,ncm, c *q1,igridsec,fn) c - Subroutine to make GMT calls to do a B/W vector plot c - Updated 10/2/2015 to allow this subroutine to work c - earlier (in makeplotfiles01.f), before "igridsec" was defined. c - See DRU-11, p. 139 implicit real*8(a-h,o-z) integer*4 maxplots character*3 ele,elecap character*200 fname character*10 olddtm,newdtm,region real*8 bw(maxplots),be(maxplots),bs(maxplots),bn(maxplots) real*4 jm(maxplots) real*4 b1(maxplots),b2(maxplots) character*10 fn(maxplots) character*5 extra character*20 units character*20 gridnote character*3 filter integer*4 igridsec,igridsec2 real*8 lorvopc,lorvog c ---------------------------------------------------- c - All vector type files begin with either "vm" c - (meters) or "vs" (arcseconds). c - The 3rd element tells me all, thinned, dropped or RMS c - The 4th/5th elements tell me if these are coordinate c - differnces or double differences. Elements 6,7,8 c - tell me what data we're plotting. c ---------------------------------------------------- c------------------------------------------------------------------- c - New for VERTCON3: c Elements 1-2: No longer have a "vs" c Elements 4-5: Can have "cd" or "dd" still, but new codes for NGVD29/NAVD88/CONUS are: c v2 = Interpolated H88-H29 from VERTCON 2 grid c rz = Diff btwen "cd" and "v2" c Elements 6-8: Type "oht" is added as the only one allowed c------------------------------------------------------------------- c ---------------------------------------------------- c - FAILSAFES: BEGIN c ---------------------------------------------------- if(fname(3:3).eq.'t')then extra='-thin' elseif(fname(3:3).eq.'d')then extra='-drop' elseif(fname(3:3).eq.'a')then extra='-all ' elseif(fname(3:3).eq.'r')then extra='-RMSd' else write(6,1)trim(fname) stop endif 1 format('FATAL in coplotvc. Bad character in spot 3: ',a) c - 2019 02 15 c if(.not.(fname(1:2).eq.'vm'.or.fname(1:2).eq.'vs'))then if(.not.(fname(1:2).eq.'vm'))then write(6,2)trim(fname) stop endif 2 format('FATAL in coplotvc. Bad character in spots 1-2: ',a) c - 2019 02 15 c if(.not.(fname(4:5).eq.'cd' .or. fname(4:5).eq.'dd' .or. c * fname(4:5).eq.'gi' ))then if(.not.(fname(4:5).eq.'cd' .or. fname(4:5).eq.'dd' .or. * fname(4:5).eq.'rz' .or. fname(4:5).eq.'v2' .or. * fname(4:5).eq.'gi' ))then write(6,3)trim(fname) stop endif 3 format('FATAL in coplotvc. Bad character in spots 4-5: ',a) if(fname(6:8).ne.ele)then write(6,4)trim(fname),ele stop endif 4 format('FATAL in coplotvc. Bad match of fname / ele: ',a,1x,a) c - Just in case we forgot to set "igridsec" to be -1 c - when we came here from makeplotfiles01: if(fname(1:3).eq.'vmacd' )then igridsec = -1 endif if (fname(2:2).eq.'m')then units='meters' elseif(fname(2:2).eq.'s')then units='arcseconds' else write(6,5)trim(fname) stop endif 5 format('FATAL in coplotvc. Bad units in name: ',a) if(igridsec.le.0)then gridnote = '' else write(gridnote,10)igridsec endif 10 format('(',i0,' sec)') c ---------------------------------------------------- c - FAILSAFES: END c ---------------------------------------------------- c ---------------------------------------------- c - GMT COMMANDS: BEGIN c ---------------------------------------------- c - Header of commands/echoes: write(99,991)ele,extra,trim(units),trim(region),trim(fn(ij)), *ele,extra,trim(units),trim(region),trim(fn(ij)) c - Make the color palette: write(99,901)cptlo,cpthi,cptin 901 format( *'gmtset FORMAT_FLOAT_OUT %.12G',/, *'makecpt -Crainbow -T',f0.10,'/',f0.10,'/',f0.10, *' > temp.cpt',/, *'gmtset FORMAT_FLOAT_OUT %.3G') c - Plot the actual vectors, and the title, with PSXY command write(99,903)trim(fname),bw(ij),be(ij),bs(ij),bn(ij), *jm(ij),b1(ij),b2(ij),trim(newdtm),trim(olddtm),elecap, *extra,trim(gridnote),trim(region),trim(fn(ij)),fname(2:5),filter 903 format('psxy ',a,' -Xc -R',f0.10,'/',f0.10,'/',sp,f0.10,'/',f0.10, * ss,' -JM',f3.1,'i -B',f0.10,'/',f0.10,':."', * 'VERTCON v3.0 ',a,' minus ',a,' ',a3,a5,a, * 1x,a,'-',a,1x,a,'-f',a3, * '": -Sc0.01i -Ctemp.cpt -K > plot.ps') c - Put in the color bar qqa = jm(ij)/2 write(99,904)qqa,trim(units) c - This version, without the "-Ba" option, will label the scale bar c - exactly as per the CPT, but the problem is that sometimes c - the "0" in the CPT is really like "10e-16", which then throws c - 16 digits on everything on the scale bar and you can't read it. c - This is a problem inside of "makecpt" and I couldn't solve it. c - As such, I have the next version of 904 which uses "-Ba" to c - automatically pick labels in psscale. c 904 format( c *'psscale -D',f3.1,'i/-0.4i/4.0i/0.2ih', c *' -Ctemp.cpt -I0.4 -By+l',a, c *': -O -K >> plot.ps') 904 format( *'psscale -D',f3.1,'i/-0.4i/4.0i/0.2ih', *' -Ctemp.cpt -I0.4 -Ba -By+l',a, *' -O -K >> plot.ps') c - ALWAYS plot the coast last, as it closes the PS file call plotcoast(region,99) c - Convert PS to JPG write(99,905) c - Rename the JPG to my naming scheme write(99,906)trim(fname),trim(fn(ij)) c ---------------------------------------------- c - GMT COMMANDS: END c ---------------------------------------------- 991 format( *'# -----------------------------------------------------',/, *'# vectors in ',a,a,1x,a,1x,a,1x,a,/, *'# -----------------------------------------------------',/, *'echo ...vectors in ',a,a,1x,a,1x,a,1x,a) c - 905 = Convert PS to JPG 905 format('psconvert plot.ps -Tj -P -A ') c - Rename to our naming scheme 906 format('mv -f plot.jpg ',a, *'.',a,'.jpg',/,'rm -f plot.ps',/, *'rm -f temp.cpt') return end