City temps on maps like Weatherbell

All other topics about postprocessing model data (GrADS and other software), about other numerical weather prediction software (including WRF-NMM and WRF-ARW discussion unrelated to UEMS/WRF EMS), and general meteorology talk go in this forum.
Post Reply
Cory
Posts: 10
Joined: Fri Sep 06, 2013 10:43 pm

City temps on maps like Weatherbell

Post by Cory » Sat Sep 28, 2013 2:23 pm

Hi everyone...I'm running the Kumar GrADS scripts and I'm wondering if anyone here has figured out how to plot the modeled values at different locations like Ryan at Weatherbell. I've been using gxout grid to plot grid temps, but having values at certain locations looks cleaner. Any help is appreciated.

Image

meteoadriatic
Posts: 1584
Joined: Wed Aug 19, 2009 10:05 am

Re: City temps on maps like Weatherbell

Post by meteoadriatic » Sat Sep 28, 2013 3:32 pm

Try this:

In the script, after you draw shaded background or whatever you do...

Code: Select all

'set strsiz 0.09'
'set string 1 c 1'
'set gxout print'

# Location xx yy :
'set lat xx' ; 'set lon yy'
'd tmp2m-273'
temp=sublin(result,2)
'query w2xy xx yy'
x = subwrd(result,3)
y = subwrd(result,6)
'draw string 'x' 'y' 'temp

 ... then put same block for each other location ...

return
Change xx to your latitude, yy to your longitude.

You can write this little bit nicer using variables instead of hardcoded lat and lon.

Cory
Posts: 10
Joined: Fri Sep 06, 2013 10:43 pm

Re: City temps on maps like Weatherbell

Post by Cory » Mon Sep 30, 2013 7:30 pm

Thanks for the response! I can't seem to get that working. I put it after gxout shaded, at various other points, and just before the image creation line, and all I get is a white map with the title and valid time.

I'm hoping it is not something to do with time loops, as the preconfigured Kumar scripts don't look very friendly for alteration.

meteoadriatic
Posts: 1584
Joined: Wed Aug 19, 2009 10:05 am

Re: City temps on maps like Weatherbell

Post by meteoadriatic » Mon Sep 30, 2013 8:37 pm

Let's see your whole script?

Cory
Posts: 10
Joined: Fri Sep 06, 2013 10:43 pm

Re: City temps on maps like Weatherbell

Post by Cory » Tue Oct 08, 2013 2:14 am

Here's the full section of one field in the script file (2m temp). It's a modified form of the Kumar scripts that come with EMS, so it exists as a section in the temp.gs script in the /util/grads/products/kumar/scripts directory.

Code: Select all

*-----------------------2M TEMP-----------------------------------------------------
t=1
fh=0
while (t<=nf)
'set t 't
cdate = mydate()
chr = subwrd(cdate,1)
cdy = subwrd(cdate,2)
cmo = subwrd(cdate,3)
cyr = subwrd(cdate,4)
*
'set vpage 0 11 1 8.5'
'set grads off'
'set gxout shaded'
'set rgb 33 248 50 60'
'set rgb 34 255 50 89'
'set rgb 35 255 50 185'
'set rgb 36 248 50 255'
'set rgb 37 224 50 255'
'set rgb 38 195 50 255'
'set rgb 39 175 50 255'
'set rgb 40 161 50 255'
'set rgb 41 137 50 255'
'set rgb 42 118 74 255'
'set rgb 43 98 74 255'
'set rgb 44 79 50 255'
'set rgb 45 50 50 255'
'set rgb 46 50 74 255'
'set rgb 47 50 89 255'
'set rgb 48 50 113 255'
'set rgb 49 50 146 255'
'set rgb 50 50 175 255'
'set rgb 51 50 204 255'
'set rgb 52 50 224 255'
'set rgb 53 50 255 253'
'set rgb 54 50 255 228'
'set rgb 55 50 255 200'
'set rgb 56 50 255 161'
'set rgb 57 50 255 132'
'set rgb 58 50 255 103'
'set rgb 59 50 255 79'
'set rgb 60 50 255 60'
'set rgb 61 79 255 50'
'set rgb 62 132 255 50'
'set rgb 63 171 255 50'
'set rgb 64 204 255 50'
'set rgb 65 224 255 50'
'set rgb 66 253 255 50'
'set rgb 67 255 233 50'
'set rgb 68 255 224 50'
'set rgb 69 255 209 50'
'set rgb 70 255 204 50'
'set rgb 71 255 185 50'
'set rgb 72 255 161 50'
'set rgb 73 255 146 50'
'set rgb 74 255 127 50'
'set rgb 75 255 118 50'
'set rgb 76 255 93 50'
'set rgb 77 255 74 50'
'set rgb 78 255 60 50'
'set rgb 79 255 33 33'
'set rgb 80 255 0 0'
'set rgb 81 235 10 0'
'set rgb 82 215 20 0'
'set rgb 83 195 30 0'
'set rgb 84 175 40 0'
'set rgb 85 165 45 0'
'set rgb 86 155 50 0'
'set rgb 87 145 55 0'
'set rgb 88 135 60 0'
'set rgb 89 120 60 0'
'set rgb 90 100 60 0'
'set rgb 91 80 60 0'
'set rgb 20 250 240 230'
'set rgb 21 240 220 210'
'set rgb 22 225 190 180'
'set rgb 23 200 160 150'
'set rgb 24 180 140 130'
'set rgb 25 160 120 110'
'set rgb 26 140 100 90'
'set clevs 0 2 4 6 8 10 12 14 16 18 20 22 24 26 28 30 32 34 36 38 40 42 44 46 48 50 52 54 56 58 60 62 64 66 68 70 72 74 76 78 80 82 84 86 88 90 92 94 96 98 100 102 104 106 108 110'
'set ccols 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 85 86 87 88 89 90 91'
'd (TMP2m-273)*9/5+32'
'set gxout contour'
'set rgb 29 245 117 23'
'set clevs 32'
'set ccols 1'
'set clab on'
'set cthick 12'
'd (TMP2m-273)*9/5+32'
*'run cbarn.gs'
'run cbarntemp'
'define temp = (TMP2m-273)*9/5+32'
'set gxout grid' 
'set gridln off'
'set ccolor 1' 
'set cthick 1'
'set digsiz 0.11' 
'd skip(temp,8,8)'
'set strsiz 0.07'
'set string 15 c 6'
'set string 1 c 1 0'
'set strsiz 0.2'
'draw shp /usr1/wrfems/util/grads/shapefiles/c_15au13.shp'
'draw string 5.5 7.3  2-meter Temperature (F)'
'set strsiz 0.13'
'draw string 5.5 6.9 Forecast Valid  'chr''cdy''cmo''cyr''
'printim  'fh'_2mtemp.png x1000 y800' 
'clear'
t=t+1
fh=fh+fint
endwhile

meteoadriatic
Posts: 1584
Joined: Wed Aug 19, 2009 10:05 am

Re: City temps on maps like Weatherbell

Post by meteoadriatic » Tue Oct 08, 2013 8:56 am

It was my mistake. In query w2xy lon should be first, lat second. Also with few other little changes:

Code: Select all

# This draws temperature at 49N, 10E:
######################################
'set strsiz 0.13'
'set string 1 c 6'
'set gxout print'
'set prnopts %0.0f'

# Location 49 10 :
'set lat 49' ; 'set lon 10'
'd (TMP2m-273)*9/5+32'
temp=sublin(result,3)
'query w2xy 10 49'
x = subwrd(result,3)
y = subwrd(result,6)
'draw string 'x' 'y' 'temp
#######################################
Whole test script based on yours:

Code: Select all

'reinit'
'set display color white'
'set mpdset hires'
'c'
'open wrf_d01.ctl'
'set grads off'
'set t 2'
'set gxout shaded'
'set rgb 33 248 50 60'
'set rgb 34 255 50 89'
'set rgb 35 255 50 185'
'set rgb 36 248 50 255'
'set rgb 37 224 50 255'
'set rgb 38 195 50 255'
'set rgb 39 175 50 255'
'set rgb 40 161 50 255'
'set rgb 41 137 50 255'
'set rgb 42 118 74 255'
'set rgb 43 98 74 255'
'set rgb 44 79 50 255'
'set rgb 45 50 50 255'
'set rgb 46 50 74 255'
'set rgb 47 50 89 255'
'set rgb 48 50 113 255'
'set rgb 49 50 146 255'
'set rgb 50 50 175 255'
'set rgb 51 50 204 255'
'set rgb 52 50 224 255'
'set rgb 53 50 255 253'
'set rgb 54 50 255 228'
'set rgb 55 50 255 200'
'set rgb 56 50 255 161'
'set rgb 57 50 255 132'
'set rgb 58 50 255 103'
'set rgb 59 50 255 79'
'set rgb 60 50 255 60'
'set rgb 61 79 255 50'
'set rgb 62 132 255 50'
'set rgb 63 171 255 50'
'set rgb 64 204 255 50'
'set rgb 65 224 255 50'
'set rgb 66 253 255 50'
'set rgb 67 255 233 50'
'set rgb 68 255 224 50'
'set rgb 69 255 209 50'
'set rgb 70 255 204 50'
'set rgb 71 255 185 50'
'set rgb 72 255 161 50'
'set rgb 73 255 146 50'
'set rgb 74 255 127 50'
'set rgb 75 255 118 50'
'set rgb 76 255 93 50'
'set rgb 77 255 74 50'
'set rgb 78 255 60 50'
'set rgb 79 255 33 33'
'set rgb 80 255 0 0'
'set rgb 81 235 10 0'
'set rgb 82 215 20 0'
'set rgb 83 195 30 0'
'set rgb 84 175 40 0'
'set rgb 85 165 45 0'
'set rgb 86 155 50 0'
'set rgb 87 145 55 0'
'set rgb 88 135 60 0'
'set rgb 89 120 60 0'
'set rgb 90 100 60 0'
'set rgb 91 80 60 0'
'set rgb 20 250 240 230'
'set rgb 21 240 220 210'
'set rgb 22 225 190 180'
'set rgb 23 200 160 150'
'set rgb 24 180 140 130'
'set rgb 25 160 120 110'
'set rgb 26 140 100 90'
'set clevs 0 2 4 6 8 10 12 14 16 18 20 22 24 26 28 30 32 34 36 38 40 42 44 46 48 50 52 54 56 58 60 62 64 66 68 70 72 74 76 78 80 82 84 86 88 90 92 94 96 98 100 102 104 106 108 110'
'set ccols 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 85 86 87 88 89 90 91'
'd (TMP2m-273)*9/5+32'
'set gxout contour'
'set rgb 29 245 117 23'
'set clevs 32'
'set ccols 1'
'set clab on'
'set cthick 6'
'd (TMP2m-273)*9/5+32'
'run cbarn.gs'
'define temp = (TMP2m-273)*9/5+32'
#'set gxout grid'
#'set gridln off'
#'set ccolor 1'
#'set cthick 1'
#'set digsiz 0.11'
#'d skip(temp,8,8)'

# This draws temperature at 49N, 10E:
######################################
'set strsiz 0.13'
'set string 1 c 6'
'set gxout print'
'set prnopts %0.0f'

# Location 49 10 :
'set lat 49' ; 'set lon 10'
'd (TMP2m-273)*9/5+32'
temp=sublin(result,3)
'query w2xy 10 49'
x = subwrd(result,3)
y = subwrd(result,6)
'draw string 'x' 'y' 'temp
#######################################

'set strsiz 0.07'
'set string 15 c 6'
'set string 1 c 1 0'
'set strsiz 0.2'
'draw string 5.5 7.3  2-meter Temperature (F)'
'set strsiz 0.13'
'draw string 5.5 6.9 Forecast Valid  'chr''cdy''cmo''cyr''
'printim  2mtemp.png x1000 y800' 
Screenshot of running it:
http://www.meteoadriatic.net/pub/wrf/weatherbell.png

I hope this helps!

BCHurricane89
Posts: 38
Joined: Thu Apr 28, 2011 4:20 am

Re: City temps on maps like Weatherbell

Post by BCHurricane89 » Fri Oct 11, 2013 11:00 pm

There has to be a way to do it without have to do the set lat/lon command as that screws up the lat/lon setting alltogether. The way I run my scripts is to have a main script that opens the files, sets the vpage, parea, labels, lat/lon for images, etc. From here I just have a listing of every script I want to run, say 2mtemp.gs. It then goes in, and does the maps. The way you have written the script, it set the lat/lon again, so when my script goes to the next timestep, nothing gets plotted because the lat/lon has been set based on the portion of the script that is suppose to 'stamp' the temperature on the maps.

Could you not simply do something like:

Code: Select all

# This draws temperature at 49N, 10E:
######################################
'set strsiz 0.13'
'set string 1 c 6'
'set gxout print'
'set prnopts %0.0f'

# Location 49 10 :
'query w2xy 10 49'
x = subwrd(result,3)
y = subwrd(result,6)
'draw string 'x' 'y' 'temp
#######################################
Assuming temp is already defined.

Cory
Posts: 10
Joined: Fri Sep 06, 2013 10:43 pm

Re: City temps on maps like Weatherbell

Post by Cory » Sun Jan 26, 2014 12:53 am

BCHurricane89 wrote:There has to be a way to do it without have to do the set lat/lon command as that screws up the lat/lon setting alltogether. The way I run my scripts is to have a main script that opens the files, sets the vpage, parea, labels, lat/lon for images, etc. From here I just have a listing of every script I want to run, say 2mtemp.gs. It then goes in, and does the maps. The way you have written the script, it set the lat/lon again, so when my script goes to the next timestep, nothing gets plotted because the lat/lon has been set based on the portion of the script that is suppose to 'stamp' the temperature on the maps.

Could you not simply do something like:

Code: Select all

# This draws temperature at 49N, 10E:
######################################
'set strsiz 0.13'
'set string 1 c 6'
'set gxout print'
'set prnopts %0.0f'

# Location 49 10 :
'query w2xy 10 49'
x = subwrd(result,3)
y = subwrd(result,6)
'draw string 'x' 'y' 'temp
#######################################
Assuming temp is already defined.
Using this method is the first time I've gotten something to plot on each graphic (the previous one left all graphics after the first one blank). I tried to define "temp" using "define temp" using various combinations, but this only plots the actual word "temp" at the defined lat/lon. I think meteo's method would have worked but my scripts use time loops and I think that's the problem. This one seems like a step in the right direction.

Cory
Posts: 10
Joined: Fri Sep 06, 2013 10:43 pm

Re: City temps on maps like Weatherbell

Post by Cory » Mon Jan 27, 2014 4:58 am

Found a solution!

Code: Select all

'set gxout print'
'set prnopts %0.0f'
lon = -99.46
lat = 29.32
'q w2xy 'lon' 'lat
say result
x = subwrd(result,3)
y = subwrd(result,6)
'define t2 = (TMP2m-273)*9/5+32'
'd t2'
'q defval t2 'x' 'y
t2m = subwrd(result,3)
t2m = math_nint(t2m)
'draw string 'x' 'y' 't2m

Cory
Posts: 10
Joined: Fri Sep 06, 2013 10:43 pm

Re: City temps on maps like Weatherbell

Post by Cory » Mon Dec 01, 2014 2:36 am

To wrap this up for others having problems: my previous solution wasn't quite a "solution", because it turned out that it was plotting temps, but had the incorrect values. I finally figured out it was because (at least in my case), it was grabbing the temp based on the w2xy values, rather than gridpoints. To fix that, I had to use two queries: one to get the lat/lon to position the value on the map, and the other to get the actual temp value at the nearest point (i,j rather than x,y). So, with this method, the correct temperature from the nearest gridpoint (i,j) is plotted on the map at x,y.

Code: Select all

***First query to get x,y from lat/lon:
lon = -99.0855
lat = 29.9767
'q w2xy 'lon' 'lat
x = subwrd(result,3)
y = subwrd(result,6)
***Second query to get temp at i,j rather than x,y: 
'q w2gr  'lon' 'lat
x2 = subwrd(result,3)
y2 = subwrd(result,6)
'define t2 = (TMP2m-273.15)*1.8+32'
***Query (defval) value of temp at i,j
'q defval t2 'x2' 'y2
t2m = subwrd(result,3)
**Round temperature value:
t2m = math_nint(t2m)
'draw string 'x' 'y' 't2m

Post Reply