1--[[ 2 plshade demo, using color fill. 3 4 Copyright (C) 2008 Werner Smekal 5 6 This file is part of PLplot. 7 8 PLplot is free software you can redistribute it and/or modify 9 it under the terms of the GNU Library General Public License as published 10 by the Free Software Foundation either version 2 of the License, or 11 (at your option) any later version. 12 13 PLplot is distributed in the hope that it will be useful, 14 but WITHOUT ANY WARRANTY without even the implied warranty of 15 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 16 GNU Library General Public License for more details. 17 18 You should have received a copy of the GNU Library General Public License 19 along with PLplot if not, write to the Free Software 20 Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA 21--]] 22 23-- initialise Lua bindings for PLplot examples. 24dofile("plplot_examples.lua") 25 26-- Fundamental settings. See notes[] for more info. 27ns = 20 -- Default number of shade levels 28nx = 35 -- Default number of data points in x 29ny = 46 -- Default number of data points in y 30exclude = 0 -- By default do not plot a page illustrating 31 -- exclusion. API is probably going to change 32 -- anyway, and cannot be reproduced by any 33 -- front end other than the C one. 34 35-- polar plot data 36PERIMETERPTS = 100 37 38-- Transformation function 39tr = {} 40 41function mypltr(x, y) 42 tx = tr[1] * x + tr[2] * y + tr[3] 43 ty = tr[4] * x + tr[5] * y + tr[6] 44 45 return tx, ty 46end 47 48---------------------------------------------------------------------------- 49-- f2mnmx 50-- 51-- Returns min & max of input 2d array. 52---------------------------------------------------------------------------- 53function f2mnmx(f, nx, ny) 54 fmax = f[1][1] 55 fmin = fmax 56 57 for i = 1, nx do 58 for j = 1, ny do 59 fmax = math.max(fmax, f[i][j]) 60 fmin = math.min(fmin, f[i][j]) 61 end 62 end 63 64 return fmin, fmax 65end 66 67 68function zdefined(x, y) 69 z = math.sqrt(x^2 + y^2) 70 71 return z<0.4 or z>0.6 72end 73 74-- return single bit (for OR) 75function bit(x,b) 76 return ((x % 2^b) - (x % 2^(b-1)) > 0) 77end 78 79-- logic OR for number values 80function lor(x,y) 81 result = 0 82 for p=1,16 do result = result + (((bit(x,p) or bit(y,p)) == true) and 2^(p-1) or 0) end 83 return result 84end 85 86 87---------------------------------------------------------------------------- 88-- main 89-- 90-- Does several shade plots using different coordinate mappings. 91---------------------------------------------------------------------------- 92 93px = {} 94py = {} 95 96fill_width = 2. 97cont_color = 0 98cont_width = 0. 99 100axis_opts = { "bcvtm" } 101num_values = {} 102values = {} 103axis_ticks = { 0.0 } 104axis_subticks = { 0 } 105label_opts = { pl.PL_COLORBAR_LABEL_BOTTOM } 106labels = { "Magnitude" } 107 108-- Parse and process command line arguments 109pl.parseopts(arg, pl.PL_PARSE_FULL) 110 111-- Load colour palettes 112pl.spal0("cmap0_black_on_white.pal"); 113pl.spal1("cmap1_gray.pal",1); 114 115-- Reduce colors in cmap 0 so that cmap 1 is useful on a 16-color display 116pl.scmap0n(3) 117 118-- Initialize plplot 119pl.init() 120 121-- Set up transformation function 122tr = { 2/(nx-1), 0, -1, 0, 2/(ny-1), -1 } 123 124-- Allocate data structures 125clevel = {} 126shedge = {} 127z = {} 128w = {} 129 130-- Set up data array 131for i = 1, nx do 132 x = (i-1 - math.floor(nx/2))/math.floor(nx/2) 133 z[i] = {} 134 w[i] = {} 135 for j = 1, ny do 136 y = (j-1 - math.floor(ny/2))/math.floor(ny/2)-1 137 z[i][j] = -math.sin(7*x) * math.cos(7*y) + x^2 - y^2 138 w[i][j] = -math.cos(7*x) * math.sin(7*y) + 2*x*y 139 end 140end 141 142zmin, zmax = f2mnmx(z, nx, ny) 143for i = 1, ns do 144 clevel[i] = zmin + (zmax-zmin)*(i-0.5)/ns 145end 146 147for i = 1, ns+1 do 148 shedge[i] = zmin + (zmax-zmin)*(i-1)/ns 149end 150 151-- Set up coordinate grids 152cgrid1 = {} 153cgrid1["xg"] = {} 154cgrid1["yg"] = {} 155cgrid1["nx"] = nx 156cgrid1["ny"] = ny 157 158cgrid2 = {} 159cgrid2["xg"] = {} 160cgrid2["yg"] = {} 161cgrid2["nx"] = nx 162cgrid2["ny"] = ny 163 164for i = 1, nx do 165 cgrid2["xg"][i] = {} 166 cgrid2["yg"][i] = {} 167 for j = 1, ny do 168 x, y = mypltr(i-1, j-1) 169 170 argx = x*math.pi/2 171 argy = y*math.pi/2 172 distort = 0.4 173 174 cgrid1["xg"][i] = x + distort * math.cos(argx) 175 cgrid1["yg"][j] = y - distort * math.cos(argy) 176 177 cgrid2["xg"][i][j] = x + distort * math.cos(argx) * math.cos(argy) 178 cgrid2["yg"][i][j] = y - distort * math.cos(argx) * math.cos(argy) 179 end 180end 181 182-- Plot using identity transform 183pl.adv(0) 184pl.vpor(0.1, 0.9, 0.1, 0.9) 185pl.wind(-1, 1, -1, 1) 186 187pl.psty(0) 188 189pl.shades(z, -1, 1, -1, 1, shedge, fill_width, cont_color, cont_width, 1) 190 191-- Smaller text 192pl.schr( 0.0, 0.75 ) 193-- Small ticks on the vertical axis 194pl.smaj( 0.0, 0.5 ) 195pl.smin( 0.0, 0.5 ) 196 197num_values[1] = ns + 1 198values[1] = shedge 199colorbar_width, colorbar_height = pl.colorbar( lor(pl.PL_COLORBAR_SHADE, pl.PL_COLORBAR_SHADE_LABEL), 0, 0.005, 0.0, 0.0375, 0.875, 0, 1, 1, 0.0, 0.0, cont_color, cont_width, label_opts, labels, axis_opts, axis_ticks, axis_subticks, num_values, values ) 200 201-- Reset text and tick sizes 202pl.schr( 0.0, 1.0 ) 203pl.smaj( 0.0, 1.0 ) 204pl.smin( 0.0, 1.0 ) 205 206pl.col0(1) 207pl.box("bcnst", 0, 0, "bcnstv", 0, 0) 208pl.col0(2) 209 210--pl.cont(w, 1, nx, 1, ny, clevel, mypltr, {}) 211pl.lab("distance", "altitude", "Bogon density") 212 213-- Plot using 1d coordinate transform 214 215-- Load colour palettes 216pl.spal0("cmap0_black_on_white.pal"); 217pl.spal1("cmap1_blue_yellow.pal",1); 218 219-- Reduce colors in cmap 0 so that cmap 1 is useful on a 16-color display 220pl.scmap0n(3); 221 222pl.adv(0) 223pl.vpor(0.1, 0.9, 0.1, 0.9) 224pl.wind(-1, 1, -1, 1) 225 226pl.psty(0) 227 228pl.shades(z, -1, 1, -1, 1, shedge, fill_width, cont_color, cont_width, 1, "pltr1", cgrid1) 229 230-- Smaller text 231pl.schr( 0.0, 0.75 ) 232-- Small ticks on the vertical axis 233pl.smaj( 0.0, 0.5 ) 234pl.smin( 0.0, 0.5 ) 235 236num_values[1] = ns + 1 237values[1] = shedge 238 239colorbar_width, colorbar_height = pl.colorbar( lor(pl.PL_COLORBAR_SHADE, pl.PL_COLORBAR_SHADE_LABEL), 0, 0.005, 0.0, 0.0375, 0.875, 0, 1, 1, 0.0, 0.0, cont_color, cont_width, label_opts, labels, axis_opts, axis_ticks, axis_subticks, num_values, values ) 240 241-- Reset text and tick sizes 242pl.schr( 0.0, 1.0 ) 243pl.smaj( 0.0, 1.0 ) 244pl.smin( 0.0, 1.0 ) 245 246pl.col0(1) 247pl.box("bcnst", 0, 0, "bcnstv", 0, 0) 248pl.col0(2) 249pl.lab("distance", "altitude", "Bogon density") 250 251-- Plot using 2d coordinate transform 252 253-- Load colour palettes 254pl.spal0("cmap0_black_on_white.pal"); 255pl.spal1("cmap1_blue_red.pal",1); 256 257-- Reduce colors in cmap 0 so that cmap 1 is useful on a 16-color display 258pl.scmap0n(3); 259 260pl.adv(0) 261pl.vpor(0.1, 0.9, 0.1, 0.9) 262pl.wind(-1, 1, -1, 1) 263 264pl.psty(0) 265 266pl.shades(z, -1, 1, -1, 1, shedge, fill_width, cont_color, cont_width, 0, "pltr2", cgrid2) 267 268-- Smaller text 269pl.schr( 0.0, 0.75 ) 270-- Small ticks on the vertical axis 271pl.smaj( 0.0, 0.5 ) 272pl.smin( 0.0, 0.5 ) 273 274num_values[1] = ns + 1 275values[1] = shedge 276colorbar_width, colorbar_height = pl.colorbar( lor(pl.PL_COLORBAR_SHADE, pl.PL_COLORBAR_SHADE_LABEL), 0, 0.005, 0.0, 0.0375, 0.875, 0, 1, 1, 0.0, 0.0, cont_color, cont_width, label_opts, labels, axis_opts, axis_ticks, axis_subticks, num_values, values ) 277 278-- Reset text and tick sizes 279pl.schr( 0.0, 1.0 ) 280pl.smaj( 0.0, 1.0 ) 281pl.smin( 0.0, 1.0 ) 282 283pl.col0(1) 284pl.box("bcnst", 0, 0, "bcnstv", 0, 0) 285pl.col0(2) 286pl.cont(w, 1, nx, 1, ny, clevel, "pltr2", cgrid2) 287 288pl.lab("distance", "altitude", "Bogon density, with streamlines") 289 290-- Plot using 2d coordinate transform 291 292-- Load colour palettes 293pl.spal0(""); 294pl.spal1("",1); 295 296-- Reduce colors in cmap 0 so that cmap 1 is useful on a 16-color display 297pl.scmap0n(3); 298 299pl.adv(0) 300pl.vpor(0.1, 0.9, 0.1, 0.9) 301pl.wind(-1, 1, -1, 1) 302 303pl.psty(0) 304 305pl.shades(z, -1, 1, -1, 1, shedge, fill_width, 2, 3., 0, "pltr2", cgrid2) 306 307-- Smaller text 308pl.schr( 0.0, 0.75 ) 309-- Small ticks on the vertical axis 310pl.smaj( 0.0, 0.5 ) 311pl.smin( 0.0, 0.5 ) 312 313num_values[1] = ns + 1 314values[1] = shedge 315colorbar_width, colorbar_height = pl.colorbar( lor(pl.PL_COLORBAR_SHADE, pl.PL_COLORBAR_SHADE_LABEL), 0, 0.005, 0.0, 0.0375, 0.875, 0, 1, 1, 0.0, 0.0, 2, 3., label_opts, labels, axis_opts, axis_ticks, axis_subticks, num_values, values ) 316 317-- Reset text and tick sizes 318pl.schr( 0.0, 1.0 ) 319pl.smaj( 0.0, 1.0 ) 320pl.smin( 0.0, 1.0 ) 321 322pl.col0(1) 323pl.box("bcnst", 0, 0, "bcnstv", 0, 0) 324pl.col0(2) 325 326pl.lab("distance", "altitude", "Bogon density") 327 328-- Note this exclusion API will probably change. 329 330-- Plot using 2d coordinate transform and exclusion 331if exclude~=0 then 332 333 -- Load colour palettes 334 pl.spal0("cmap0_black_on_white.pal"); 335 pl.spal1("cmap1_gray.pal",1); 336 337 -- Reduce colors in cmap 0 so that cmap 1 is useful on a 16-color display 338 pl.scmap0n(3); 339 340 pl.adv(0) 341 pl.vpor(0.1, 0.9, 0.1, 0.9) 342 pl.wind(-1, 1, -1, 1) 343 344 plpsty(0) 345 346 pl.shades(z, zdefined, -1, 1, -1, 1, shedge, fill_width, cont_color, cont_width, 347 0, "pltr2", cgrid2) 348 349 pl.col0(1) 350 pl.box("bcnst", 0, 0, "bcnstv", 0, 0) 351 352 pl.lab("distance", "altitude", "Bogon density with exclusion") 353end 354 355-- Example with polar coordinates. 356 357-- Load colour palettes 358pl.spal0("cmap0_black_on_white.pal"); 359pl.spal1("cmap1_gray.pal",1); 360 361-- Reduce colors in cmap 0 so that cmap 1 is useful on a 16-color display 362pl.scmap0n(3); 363 364pl.adv(0) 365pl.vpor(.1, .9, .1, .9) 366pl.wind(-1, 1, -1, 1) 367 368pl.psty(0) 369 370-- Build new coordinate matrices. 371for i = 1, nx do 372 r = (i-1)/(nx-1) 373 for j = 1, ny do 374 t = 2*math.pi/(ny-1)*(j-1) 375 cgrid2["xg"][i][j] = r*math.cos(t) 376 cgrid2["yg"][i][j] = r*math.sin(t) 377 z[i][j] = math.exp(-r^2)*math.cos(5*math.pi*r)*math.cos(5*t) 378 end 379end 380 381-- Need a new shedge to go along with the new data set. 382zmin, zmax = f2mnmx(z, nx, ny) 383 384for i = 1, ns+1 do 385 shedge[i] = zmin + (zmax-zmin)*(i-1)/ns 386end 387 388-- Now we can shade the interior region. 389pl.shades(z, -1, 1, -1, 1, shedge, fill_width, cont_color, cont_width, 0, "pltr2", cgrid2) 390 391-- Smaller text 392pl.schr( 0.0, 0.75 ) 393-- Small ticks on the vertical axis 394pl.smaj( 0.0, 0.5 ) 395pl.smin( 0.0, 0.5 ) 396 397num_values[1] = ns + 1 398values[1] = shedge 399colorbar_width, colorbar_height = pl.colorbar( lor(pl.PL_COLORBAR_SHADE, pl.PL_COLORBAR_SHADE_LABEL), 0, 0.005, 0.0, 0.0375, 0.875, 0, 1, 1, 0.0, 0.0, cont_color, cont_width, label_opts, labels, axis_opts, axis_ticks, axis_subticks, num_values, values ) 400 401-- Reset text and tick sizes 402pl.schr( 0.0, 1.0 ) 403pl.smaj( 0.0, 1.0 ) 404pl.smin( 0.0, 1.0 ) 405 406-- Now we can draw the perimeter. (If do before, shade stuff may overlap.) 407for i = 1, PERIMETERPTS do 408 t = 2*math.pi/(PERIMETERPTS-1)*(i-1) 409 px[i] = math.cos(t) 410 py[i] = math.sin(t) 411end 412pl.col0(1) 413pl.line(px, py) 414 415-- And label the plot. 416pl.col0(2) 417pl.lab( "", "", "Tokamak Bogon Instability" ) 418 419pl.plend() 420