This page was created by the IDL library routine
mk_html_help. For more information on
this routine, refer to the IDL Online Help Navigator
or type:
? mk_html_help
at the IDL command line prompt.
Last modified: Tue Nov 06 13:52:53 2012.
NAME:
COLORBAR
PURPOSE:
The purpose of this routine is to add a color bar to the current
graphics window.
AUTHOR:
FANNING SOFTWARE CONSULTING
David Fanning, Ph.D.
1645 Sheely Drive
Fort Collins, CO 80526 USA
Phone: 970-221-0438
E-mail: davidf@dfanning.com
Coyote's Guide to IDL Programming: http://www.dfanning.com/
CATEGORY:
Graphics, Widgets.
CALLING SEQUENCE:
COLORBAR
INPUTS:
None.
KEYWORD PARAMETERS:
BOTTOM: The lowest color index of the colors to be loaded in
the bar.
CHARSIZE: The character size of the color bar annotations. Default is 1.0.
COLOR: The color index of the bar outline and characters. Default
is !P.Color..
DIVISIONS: The number of divisions to divide the bar into. There will
be (divisions + 1) annotations. The default is 6.
FONT: Sets the font of the annotation. Hershey: -1, Hardware:0, True-Type: 1.
FORMAT: The format of the bar annotations. Default is '(I5)'.
INVERTCOLORS: Setting this keyword inverts the colors in the color bar.
MAXRANGE: The maximum data value for the bar annotation. Default is
NCOLORS.
MINRANGE: The minimum data value for the bar annotation. Default is 0.
MINOR: The number of minor tick divisions. Default is 2.
NCOLORS: This is the number of colors in the color bar.
POSITION: A four-element array of normalized coordinates in the same
form as the POSITION keyword on a plot. Default is
[0.88, 0.15, 0.95, 0.95] for a vertical bar and
[0.15, 0.88, 0.95, 0.95] for a horizontal bar.
RANGE: A two-element vector of the form [min, max]. Provides an
alternative way of setting the MINRANGE and MAXRANGE keywords.
RIGHT: This puts the labels on the right-hand side of a vertical
color bar. It applies only to vertical color bars.
TITLE: This is title for the color bar. The default is to have
no title.
TOP: This puts the labels on top of the bar rather than under it.
The keyword only applies if a horizontal color bar is rendered.
VERTICAL: Setting this keyword give a vertical color bar. The default
is a horizontal color bar.
COMMON BLOCKS:
None.
SIDE EFFECTS:
Color bar is drawn in the current graphics window.
RESTRICTIONS:
The number of colors available on the display device (not the
PostScript device) is used unless the NCOLORS keyword is used.
EXAMPLE:
To display a horizontal color bar above a contour plot, type:
LOADCT, 5, NCOLORS=100
CONTOUR, DIST(31,41), POSITION=[0.15, 0.15, 0.95, 0.75], $
C_COLORS=INDGEN(25)*4, NLEVELS=25
COLORBAR, NCOLORS=100, POSITION=[0.15, 0.85, 0.95, 0.90]
MODIFICATION HISTORY:
Original version by David Fanning
Last modified by A. Mignone, 1st May 2011
(See colorbar.pro)
NAME: curl
AUTHOR: Andrea Mignone (mignone@to.astro.it)
PURPOSE: compute the curl of a vector field
SYNTAX: array = curl(u,v,x1,x2,dx1,dx2,geo=geo)
DESCRIPTION: curl returns an array of the same size
of u (or v) containing the curl
of the vector field (u,v) computed as:
d(u)/dy - d(v)/dx
The input parameters are:
u = first vector component (2-D array)
v = second vector component (2-D array)
x1 = 1-D array of n1 points containing the x1 zone centers
x2 = 1-D array of n2 points containing the x2 zone centers
dx1 = 1-D array of n1 points containing mesh widths
dx2 = 1-D array of n2 points containing mesh widths
geo = specify in which geometry the divergence
should be computed. Options are:
'cartesian', 'spherical', 'cylindrical', 'polar'
EXAMPLES
#1 Compute and display the curl of the initial magnetic field (toroidal
current) in cartesian coordinates:
IDL> display, curl(b1(0), b2(0), x1, x2, dx1, dx2, geo="cartesian"), /vbar
HISTORY: Sept 22, 2003
(See curl.pro)
NAME: DISPLAY
AUTHOR: Andrea Mignone (mignone@ph.unito.it)
PURPOSE: Make a tvscl image with axis, colorbar, and all the
nice stuffs.
SYNTAX: display, img[,x1=x1][,x2=x2][,BACKGROUND=value]
[,CBDIV=integer][,CHARSIZE=value][,COLOR=value][,/EPS]
[,FILENAME=string][,IMAX=value][,IMIN=value]
[,IMSIZE=value][,LABEL1=string][,LABEL2=string][,/NOCLOSE]
[,/SMOOTH][,TITLE=string][,/HBAR][,/VBAR][,/XSYM][,/YSYM]
[,XRANGE=[min,max]][,YRANGE=[min,max]][,NWIN=integer]
[,/POLAR][,{LFT,RGT,TOP,BOT}=value]
KEYWORDS:
x1 = a one 1-D array with the horizontal coordinates
dim(x1) = first dimension of img;
x2 = a one 1-D array with the vertical coordinates
dim(x2) = second dimension of img; both
x1 and x2 must be present one wish to
scale axis to coordinates
BACKGROUND = in integer value in the range [1,255] giving the
background color.
CBDIV = an integer number specifying the number of
divisions to divide the colorbar into.
CHARSIZE = size of characters.
COLOR = an integer in the range [1,255] giving contour and
colorbar should be > 0 and <= 255
/EPS = when this keyword is set, it produces eps output.
the image keeps the aspect ratio with a default
horizontal extent of 10cm. Use imsize to shrink/enlarge it.
FILENAME = a string giving the name of the file where output is saved;
only works when /eps is used.
IMAX = maximum floating value to which img should be scaled to;
IMIN = minimum floating value to which img should be scaled to;
IMSIZE = a magnification/shrinking scale factor;
for example, imsize=2.0 produces an image twice
as big as the original.
When graphics is NOT eps, it can also be a two
element arrays giving the new image sizes
LABEL1 = a string label for the x1 - axis
LABEL2 = a string label for the x2 - axis
NWIN = an integer selecting the window.
/POLAR = maps a polar (r,phi) domain into cartesian
coordinates. The user must also supply the
x1 (= r) and x2 (= phi) coordinates. It also
recomended that imsize = [nx,ny] be given, where
nx and ny are the final number of points in the
cartesian map.
/SMOOTH = smooth the image using cubic interpolation when resizing
TITLE = string title of the plot
/HBAR = add horizontal color bar at the bottom of the plot
/VBAR = add a vertical color bar to the right of the plot
/XSYM = symmetrize the image with respect to the x-axis
/YSYM = symmetrize the image with respect to the y-axis
XRANGE = a two-element vector specifying the abscissa of the
lower and upper boundary for the displayed image.
YRANGE = a two-element vector specifying the ordinata of the
lower and upper boundary for the displayed image.
/NOCLOSE = do not close output device (for EPS output)
LFT,RGT = these parameters can be used to specify the margin (in pixels)
BOT,TOP between the plot borders and the left (lft), right (rgt),
top (top) and bottom (bot) side of the window.
EXAMPLES:
* Example #1: Display the logarithm of the intial density and put axis,
horizontal colorbar and title:
IDL> display, alog(rho(0)), x1 = x1, x2 = x2, title = "My title", /hbar
* Example #2: Display another image in a different window with a
vertical colorbar using the red - blue color table:
IDL> loadct,33
IDL> display, img2, x1 = x1, x2 = x2, /vbar, nwin = 1
* Example #3: Display a subdomain of the original density map and leave
extra space to the right of the plot for extra annotations:
IDL> display, alog(rho(0)), x1 = x1, x2 = x2, $
xrange=[0.3,0.7], yrange=[0.4,0.8],rgt = 120
* Example #4: Mapping an image from polar to Cartesian on a 400x400 grid:
IDL> display, img, x1 = x1, x2 = x2, /polar, ims = [400,400]
* Example #5: Mapping an image from spherical to Cartesian:
IDL> display, img, x1 = x1, x2 = 0.5*!PI - x2, /polar, ims = [400,400]
LAST MODIFIED: Aug 24, 2011 by A.Mignone (mignone@ph.unito.it)
(See display.pro)
NAME: div
AUTHOR: A. Mignone
PURPOSE: to compute the divergence of a vector field
SYNTAX: array = div(u,v,x1,x2,dx1,dx2,/geo)
DESCRIPTION: div returns an array of the same size
of u (or v) containing the divergence
of the vector field (u,v) computed as:
d(Au)/dV + d(Av)/dV
The input parameters are:
u = first vector component (2-D array)
v = second vector component (2-D array)
x1 = 1-D array of n1 points containing the x1 zone centers
x2 = 1-D array of n2 points containing the x2 zone centers
dx1 = 1-D array of n1 points containing mesh widths
dx2 = 1-D array of n2 points containing mesh widths
/geo = specify in which geometry the divergence
should be computed. Options are:
/cartesian, /spherical, /cylindrical, /polar
EXAMPLES
#1 Compute and display the divergence of the magnetic field in cartesian
coordinates of the 1st output file:
IDL> display, div(b1(1), b2(1), x1, x2, dx1, dx2, /cartesian), /vbar
#2 Compute and display the initial divergence of the velocity in
cylindrical coordinates:
IDL> display, div(v1(0), v2(0), x1, x2, dx1, dx2, /cylindrical), /vbar;
HISTORY: Sept 22, 2003
(See div.pro)
NAME: div3d
AUTHOR: P. Tzeferacos, A. Mignone
PURPOSE: to compute the divergence of a vector field
SYNTAX: array = div(u,v,w,x1,x2,x3,dx1,dx2,dx3,/geo)
DESCRIPTION: The input parameters are:
u = first vector component (3-D array)
v = second vector component (3-D array)
w = third vector component (3-D array)
x1 = 1-D array of n1 points containing the x1 zone centers
x2 = 1-D array of n2 points containing the x2 zone centers
x3 = 1-D array of n3 points containing the x3 zone centers
dx1 = 1-D array of n1 points containing mesh widths
dx2 = 1-D array of n2 points containing mesh widths
dx3 = 1-D array of n3 points containing mesh widths
/geo = specify in which geometry the divergence
should be computed. Options are:
/cartesian, /spherical, /polar
HISTORY: July 12, 2007
(See div3d.pro)
NAME: extrema
AUTHOR: Andrea Mignone
DATE: April 21, 2004
PURPOSE: Find maxima/minima in a function f.
f is a 1-D vector.
SYNOPSIS: Result = extrema(f). On output
result is an integer array whose
elements are the indexes k of the
extrema of f.
KEYWORD: none
(See extrema.pro)
UNCTION FIELD_INTERP, U, V, W, x, y, z, coord
PURPOSE: Linearly interpolate the vector
fields U and V at point xp and yp
ARGUMENTS:
U,V: 2-D vectors to be interpolated
nx,ny: number of points in the two directions
x,y: 1-D abscissas and ordinates
xp,yp: location where U and V should be interpolated
EXAMPLES
#1 Plot the logarithm of the density at the initial time and overplot
2 magnetic fieldlines rooted at the (x,y) = (1.0,0.1) and (2.0,0.1)
IDL> display, alog(rho(0)), x1 = x1, x2 = x2, /vbar
IDL> field_line, b1(0), b2(0), x1, x2, 1.0, 0.1, qx, qy
IDL> oplot, qx, qy
IDL> field_line, b1(0), b2(0), x1, x2, 2.0, 0.1, qx, qy
IDL> oplot, qx, qy
(See field_line.pro)
NAME: FIELD_LINE
AUTHOR: A. Mignone (mignone@ph.unito.it)
PURPOSE: Given a 2-D vector field (U,V), it computes the field line
passing through x0, y0. Integration stops when either
the domain boundaries are reached or the max number of iteration
is exceeded.
SYNTAX: FIELD_LINE, U, V, x, y, x0, y0, xfield, yfield[,step=step]
[,method=method][,maxstep=maxstep][,minstep=minstep][,tol=tol]
ARGUMENTS:
U,V: 2-D arrays giving the 1st (U) and 2nd (V) component of the
vector field
x,y: 1-D arrays giving abscissas and ordinates
x0,y0: two scalars giving the point coordinates through which the
field line goes.
xfield,yfield: 1-D arrays containing the output abscissas and ordinates
of the field line.
KEYWORDS:
step: a scalar giving the initial step size. Default is (mean) grid spacing.
method: a string giving the integration method. The possible choices are:
"RK2" explicit, fixed step, 2nd order Runge-Kutta methods.
"RK4" explicit, fixed step, 4th order Runge-Kutta methods.
"BS23" explicit, adaptive stepsize Runge-Kutta-Fehlberg of order
3 with local truncation error based on a 2nd-order accurate
embedded solution.
"CK45" explicit, adaptive stepsize Cask-Karp of order
5 with local truncation error based on a 4th-order accurate
embedded solution.
The default is "RK2". Use an adaptive stepsize integrator
when the field line does not close on itself inside the domain.
maxstep: a scalar giving the maximum allowed integration step.
Default is 100*step.
minstep: a scalar giving the minimum allowed integration step.
Default is 0.05*step.
tol: a scalar value giving the relative tolerance during adaptive stepsize
integration. It is ignored for fixed step-size integration (such as RK2, RK4)
EXAMPLE:
* compute a field line tangent to the vector field (B1,B2) in the point
with coordinate (-1,2) using the Bogacki-Shampine method with relative
tolerance 1.e-4:
IDL> field_line, b1, b2, x1, x2, -1, 2, xl, yl, method="BS23", tol=1.e-4
IDL> oplot,xl,yl ; overplot on current window
LAST MODIFIED
September 23, 2010 by A. Mignone (mignone@ph.unito.it)
(See field_line.pro)
NAME: fourier
AUTHOR: Andrea Mignone
DATE: Feb 10, 2005
PURPOSE: Take the fourier transform of some signal and make
a log-log plot of its power in frequency space.
SYNOPSIS: fourier, time, sp, nbeg, nend[, title = title]
ARGUMENTS / KEYWORDS:
time: a 1-D vector containing the
times at which the signal
is sampled.
sp : a 1-D vector containing the
sampled data.
nbeg: the initial index of the time
window to be transformed.
In other words, only the interval
time(nbeg) < t < time(nend)
is considered.
nend: the final index.
title: an optional keyword specifying the
title of the plot
VERSION: beta
(See fourier.pro)
NAME: get_frame
AUTHOR: Andrea Mignone
PURPOSE: grab window's content and write
an image file. The output file name is required.
SYNOPSIS: get_frame,name=name[,/jpg][,/tiff][,/png][,/ppm]
[,/gif][,quality=quality]
KEYWORDS: describe the image format.
EXAMPLES
#1 Read the current window content and dump it to rho.jpg:
IDL> get_frame, name = "rho", /jpg
LAST MODIFIED: Jan 13 2006 by A. Mignone (mignone@to.astro.it)
(See get_frame.pro)
NAME: grad
AUTHOR: A. Mignone
PURPOSE: to compute the gradient of a scalar
SYNTAX: array = grad(u,v,x1,x2,dx1,dx2,/geo)
DESCRIPTION:
phi = scalar function (2-D array)
x1 = 1-D array of n1 points containing the x1 zone centers
x2 = 1-D array of n2 points containing the x2 zone centers
dx1 = 1-D array of n1 points containing mesh widths
dx2 = 1-D array of n2 points containing mesh widths
/geo = specify in which geometry the divergence
should be computed. Options are:
/cartesian, /spherical, /cylindrical, /polar
HISTORY: Sept 22, 2003
(See grad.pro)
NAME: H5LOAD
AUTHOR: C. Zanni
REQUIRES: HDF5 support.
PURPOSE: Read PLUTO variables stored in a (pixie) HDF5 data file
(static grid). It is automatically called by PLOAD when the
/H5 keyword is supplied.
SYNTAX: H5LOAD, filename
ARGUMENTS
filename = the data file name
KEYWORDS
/SILENT Set this keyword to suppress output
LAST MODIFIED
Sept 25, 2012 by C. Zanni (zanni@oato.inaf.it)
(See h5load.pro)
NAME: HDF5LOAD
AUTHOR: A. Mignone, O. Tesileanu, C. Zanni
REQUIRES: HDF5 support.
PURPOSE: Read a HDF5 Chombo data-file and store its content on the
variables vars. It is automatically called by PLOAD when the
/HDF5 keyword is supplied.
SYNTAX: HDF5LOAD, nout[, VAR=string][, FILENAME=string][, /INTERP]
[,LEVEL=value][,{X | Y | Z}RANGE=[min,max]]
[, /NO_VAR_MATCH][,/NODATA]
KEYWORDS
VAR = the variable name, e.g. "rho" or "bx2"
FILENAME = the file name. If omitted, "data.nnnn.hdf5" is assumed.
LEVEL = the output level resolution (default = 0)
INTERP = set this keyword if you wish to interpolate a coarse level
down to the resolution of the selected level using bilinear
interpolation provided by REBIN. Otherwise nearest REBIN uses
nearest neighbor sampling for both magnification and minification.
X1RANGE = a 2 element vector giving the x-range of the portion
to be extracted.
X2RANGE = a 2 element vector giving the y-range of the portion
to be extracted.
X3RANGE = a 2 element vector giving the z-range of the portion
to be extracted.
NO_VAR_MATCH = do not try to match variable names with those provides
by PLUTO and copy the corresponding storage areas
(e.g. rho, vx1, vx2, etc...). In other words, all variables
read from the HDF file are stores inside vars(i,j,k,nv).
Useful when HDF5LOAD is called by another function (e.g.
HDF5CUT).
EXAMPLES:
* Example #1: read "data.0002.hdf5" and set the output file resolution
equivalent to the 4th level of refinement:
IDL> HDF5LOAD, 2, lev=4
* Example #2: read "data.0005.hdf5" and set the output file resolution
equivalent to the 2th level of refinement:
IDL> HDF5LOAD, 5, lev=2, var="v1"
LAST MODIFIED
Oct 19, 2012 by A. Mignone (mignone@ph.unito.it) & C. Zanni (zanni@oato.inaf.it)
(See hdf5load.pro)
NAME: MIRROR
AUTHOR: A. Mignone
PURPOSE: Build a symmetric image out of a bidimensional array a.
SYNTAX: MIRROR, a[, LEFT=left][, RIGHT=right][, TOP=top][, BOTTOM=bottom]
ARGUMENTS:
a = 2D array
KEYWORDS
One keyword between LEFT, ..., BOTTOM should be given with values equal to
1 = reflect the image w/ respect to the axis
[useful w/ reflective boundary conditions]
2 = reflect the image w/ respect to the axis midpoint
[useful w/ point-reflective boundary conditions]
3 = do not reflect, simply copy the array
[useful w/ periodic boundary conditions]
LAST MODIFIED: Aug 21, 2012 by A. Mignone (mignone@ph.unito.it)
(See mirror.pro)
NAME: OPLOTBOX
AUTHOR: A. Mignone
PURPOSE: overplot the AMR box layout on the current window.
The axis of the plot should be scaled to physical units.
(e.g. display should have been called with x1=x1,x2=x2 keywords).
SYNTAX: OPLOTBOX, [, CTAB=integer][, CVAL=vector][, LRANGE=[min,max]]
[,_EXTRA=extra]
KEYWORDS
CTAB = an integer giving the color table to be used when plotting
the boxes. Default is the current color table.
When CTAB=12 or CTAB=23, color indices are pre-defined to
default values so that level 0 is always red, level 1 is green,
and so forth.
CVAL = on output, it returns an array of integers with the
corresponding color indices used to draw the boxes.
Useful for legends.
LRANGE = a two element integer array of the form [minlev, maxlev] giving
the minimum and maximum levels to be plotted.
OPLOTBOX also accepts keyword inheritance via the _EXTRA keyword to be passed to plots
EXAMPLES:
* Example #1: plot boxes on the current graphic device:
IDL> OPLOTBOX
* Example #2: display density maps leaving extra space for a
legend plot (rgt=140), overplot the boxes and then
add a legend:
IDL> display,q,x1=x1,x2=x2,rgt = 140
IDL> oplot,cval=vals
IDL> legend, ["Level "+strcompress(string(indgen(4)),/remove_all)],$
line=0,colors=vals,pos=[0.83,0.9],/normal,spac=0.1
LAST MODIFIED
May 17th, 2011 by A. Mignone (mignone@ph.unito.it)
(See oplotbox.pro)
NAME: PLOAD
AUTHOR: Andrea Mignone
PURPOSE:
Provides a convenient and flexible way to read PLUTO data files written
in different formats: binary (.dbl/.flt), hdf5 (.dbl.h5/.flt.h5) or
chombo-hdf5 (.hdf5).
Data is read and stored in the three common blocks PLUTO_GRID (grid-related
information), PLUTO_VAR (variable arrays) and PLUTO_RUN (time-step info).
PLOAD must be invoked prior to any other function for properly initializing
the common blocks containing useful information specific to the data being
read.
The common blocks are:
PLUTO_GRID -> contains grid information: this is read using
the READ_PLUTO_GRID function.
NX1 = number of points in the x1 direction
NX2 = number of points in the x2 direction
NX3 = number of points in the x3 direction
x1 = 1-D array of NX1 points containing the x1 zone centers
x2 = 1-D array of NX2 points containing the x2 zone centers
x3 = 1-D array of NX3 points containing the x3 zone centers
dx1 = 1-D array of NX1 points containing mesh widths
dx2 = 1-D array of NX2 points containing mesh widths
dx3 = 1-D array of NX3 points containing mesh widths
xpos = a 2D array containing the x-coordinate values (useful for
non-Cartesian geometry);
ypos = same as xpos but containing y-coordinate values.
AMRLevel = an array of structures, each containing a collection of
rectangular AMR boxes (see HDF5LOAD)
AMRBoxes = an array giving the level number as a function of the
coordinate indices.
PLUTO_VAR -> contains the number of variables (=NVAR) and the corresponding
file names:
rho = density
vx1, vx2, vx3 = velocity components
bx1, bx2, bx3 = magnetic field components
prs = pressure
tr1, ... = passive scalars
bx1s, bx2s, bx3s = staggered components of magnetic field (only with
double precision binary data files) available with
constrained transport.
PLUTO_RUN -> contains time-related information:
t = time at the corresponding output number
dt = time step at the corresponding output number
nlast = the last written file
File and variable names are read from the corresponding .out file (e.g.
flt.out, dbl.out, dbl.h5.out, etc...)
SYNTAX:
PLOAD, nout[, /ASSOC,][,DIR=string][, /FLOAT][,/H5][,/HDF5][,/NODATA]
[,SHRINK=integer][,/SILENT][, VAR=string][,/XYASSOC]
[{X1 | X2 | X3}RANGE=[min,max]][, _EXTRA = extra]
ARGUMENTS:
nout = the sequential output number, e.g., 0,1,2, ...
KEYWORDS:
/ASSOC = use IDL file-association rather than storing files into memory,
This can be particularly useful for large data sets.
DIR = the directory name where output files are located (e.g.
DIR = 'Data1/'). Default is current directory.
/FLOAT = set this keyword if you wish to read single precision
data rather than double.
The default will try to read dbl (if possible) and
will automatically switch to /float if dbl.out is not found.
/H5 = assume data is in (pixie) HDF5 format used by the static grid
version of PLUTO. Double precision data files (.dbl.h5) are
searched for by default.
Single precision files (.flt.h5) can be read by specifying
also the /FLOAT keyword.
/HDF5 = assume data is in a Chombo HDF5 data file.
When this keyword is specified, PLOAD only acts as a wrapper
for the HDF5LOAD procedure which actually reads the HDF5
datafile(s).
In this case, different KEYWORDs are supported, see the
documentation for HDF5LOAD.
/NODATA = read grid information only. Do not attempt to read/store
data (useful for large data sets)
SHRINK = an integer number > 1 specifying how much the data has to be
shrinked in all directions (only for .dbl and .flt data formats).
The original array is never stored into memory and resampling is
done slice by slice.
Beware that the grid is redefined accordingly.
Useful for large 3-D data sets.
IMPORTANT: non-uniform mesh are not properly accounted for
and the interpolation may be inaccurate.
/SILENT = set this keyword to suppress output generated by pload.
VAR = the name of a particular variable to be loaded, e.g., "pr" (only
for .dbl or .flt data formats). The remainig ones are ignored.
/XYASSOC = associate an array of 2-D data out of a 3-D data file (only for
.dbl and .flt data formats). Each index 'k' of the array
corresponds to an XY plane sliced at different z=z(k).
X1RANGE = a 2-element vector in the form [xbeg,xend]
giving the x1 axis range to be extracted from the data
(only for .dbl and .flt data formats).
Only the specified subset is stored into memory and the grid
is redefined accordingly. Useful for large datasets.
X2RANGE = same as x1range but for the x2-axis.
X3RANGE = same as x1range but for the x3-axis.
Please note that '/XYASSOC' and 'SHRINK' are mutually incompatible.
EXAMPLES:
* Example #1: retrieves the 10th output data set in double precision and
store variable contents into {rho, vx1, vx2, ...}:
IDL> pload, 10, DIR = 'mydata/'
* Example #2: open all files and associate the variables {rho,v1,...}
to of 2D XY slices at constant z planes.
For instance, rho(k) will be a 2D slice containing density
in the XY plane at z=z(k):
IDL> pload, /FLOAT, /XYASSOC, 7
* Example #3: open the single-precision datafile containing 'rho'
and resize the whole array by making it a factor of 2 smaller.
The new grid will be halved respect to the original one.
Other variables will be discarded:
IDL> pload, /FLOAT, 84, VAR="rho", SHRINK=2
* Example #4: extract the 3D portion of the domain specified by the range
keywords. The content is stored into {rho,vx1,vx2,...}.
IDL> pload, 1, X1RANGE=[0, 40], X2RANGE=[-10,10], X3RANGE=[-10,10]
* Example #5: load HDF5 data file "data.0004.hdf5" and store only "vx2":
IDL> pload,/hdf5,4,level=3, var="vx2"
LAST MODIFIED
Sept 25, 2012 by A. Mignone (mignone@ph.unito.it)
(See pload.pro)
NAME: polar
AUTHOR: Andrea Mignone
PURPOSE: Interpolates the surface array a from polar coordinates (r, theta) to
rectangular coordinates (x, y).
On output, a is replaced by the array in cartesian coordinates,
while r becomes x, theta becomes y. The array a can be
irregularly gridded. This function uses the polar_surface
routine of IDL.
Note: if spherical polar coordinates
are use, then simply let
x2 = !pi/2. - x2
SYNTAX: polar, a, x1, x2[,sample=sample][,missing=missing]
KEYWORD & PARAMETERS:
see IDL manual page for polar_surface.
DATE: Aug 20, 2004
(See polar.pro)
NAME: put_eps
AUTHOR: Andrea Mignone
PURPOSE: Put a single encapsulated postscript figure in
a multi-plot file.
SYNTAX: put_eps, data, x, y, pos
[,xtitle=xtitle][,ytitle=ytitle][,title=title]
[,dmax = dmax][,dmin = dmin][,/vbar][,/hbar][,color=color]
[,cbwidth=cbwidth][,charsize=charsize][,cbcharsize=cbcharsize]
[,/noxticks][,/noyticks][,xrange=xrange][,yrange=yrange]
where data is a two-dimensional array, x and y
are the coordinates, pos=[x0,y0,x1,y1] is the position
of the figure in the multiplot file,
calculated with set_multi_plot_pos.
KEYWORDS:
xtitle = the string name for the x-axis
ytitle = the string name for the y-axis
title = the string name for the title
dmax = a floating point number giving the maximum value
data value to which the image should be scaled to.
dmin = a floating point number giving the minimum value
data value to which the image should be scaled to.
vbar = set this keyword to put a vertical colorbar to the right
of the image.
hbar = set this keyword to put a horizontal colorbar to
the top of the image.
color = an integer giving the color (in the current color table)
of the axis.
charsize = size of characters used for colorbar
noxticks = suppress x-tick marks
noyticks = suppress y-tick marks
sample = increase resolution by combining sub-sampling and
interpolation.
xrange = a two element vector giving the range in the x direction
yrange = a two element vector giving the range in the y direction
COLORBAR CONTROL KEYWORDS:
cbwidth = the horizontal extent of the colorbar in units
of the image (only if /vbar or /hbar is given).
cbcharsize = size of characters used for colorbar
cbdiv = the number of divisions of the colorbar.
LAST MODIFIED: Oct 6 2010 by A. Mignone (mignone@ph.unito.it)
(See put_eps.pro)
NAME: REGRID
AUTHOR: Andrea Mignone (mignone@ph.unito.it)
PURPOSE: Linearly interpolate irregularly-gridded data to a
a regular grid. Interpolation is done
separately row by row and column by column.
On output the original arrays are replaced with the
the new ones.
SYNTAX: REGRID, img, x1, x2[,n1=integer][,n2=integer]
[, {X | Y}RANGE = [min,max]]
img = a 2-D array.
x1 = a 1-D array with the abscissa of img
x2 = a 1-D array with the ordinata of img
NOTICE: On output, img, x1 and x2 are replaced by the new
regridded arrays. x1 and x2 will have uniform spacing.
KEYWORDS:
n1 = the number of points in the regular grid (x direction)
n2 = the number of points in the regular grid (y direction)
XRANGE = a two dimensional vector [xbeg, xend] defining
the lower and upper points of the new grid in
the x direction
YRANGE = a two dimensional vector [ybeg, yend] defining
the lower and upper points of the new grid in
the y direction
LAST MODIFIED: Feb 4, 2006
(See regrid.pro)
NAME: set_multi_plot_pos
AUTHOR: Andrea Mignone (mignone@to.astro.it)
PURPOSE: Compute position coordinates for multi plot figures
SYNTAX: pos = set_multi_plot_pos(ncol,nrow),
[,top=top][,bot=bot][,lft=lft][,rgt=rgt]
[,xsep=xsep][,ysep=ysep]
where nrow and ncol are the number of rows and columns
in the plot. Note: this function only computes
the coordinates for the position vector; it does not
produce any plot!
On output, pos(i,j,*) will be a 1-D array giving
the position of the plot as (x,y,x+xsep,y+ysep).
For example, in a 2 x 2 multiplot figure you have
pos(0,0) pos(1,0)
pos(0,1) pos(1,1)
KEYWORDS:
top = a scalar specifying the space between the top margin and
the first row
bot = a scalar specifying the space between the bottom margin and
the last row
lft = a scalar specifying the space between the left margin and
the first column
rgt = a scalar specifying the space between the right margin and
the last column
xsep = horizontal separation between plots
ysep = vertical separation between plots
LAST MODIFIED: March 12, 2007
(See set_multi_plot_pos.pro)
NAME: shockfind
AUTHOR: Andrea Mignone (mignone@to.astro.it)
PURPOSE: find shocks in a multidimensional plane
SYNTAX: array = shockfind(p, vx,vy[,eps_min=eps_min][,eps_max=eps_max])
DESCRIPTION: shockfind returns a multi-D array (with dimensionality
equal to the input arrays) whose elements are
numbers between 0 (no shock) and 1 (strong shock).
A shock is detected when div.V < 0
and |grad p|/p > eps. Divergence and gradients
are undivided differences.
The strength of the shock is controlled by
eps_min < |grad p|/p < eps_max
ARGUMENTS and KEYWORDS:
p: an array giving the pressure distribution
vx: an array giving the x-component of velocity
vy: an array giving the y-component of velocity
eps_min: the minimum threshold for a shock to exist
eps_max: the maximum strength of a shock (in units
of |grad p|/p).
LAST MODIFIED: Aug 6 2006 by A.Mignone (mignone@to.astro.it)
(See shockfind.pro)
NAME: VECFIELD
AUTHOR: A. Mignone (mignone@ph.unito.it)
PURPOSE: Produces a two-dimensional velocity field plot.
SYNTAX: VECFIELD, U, V, X, Y[,SPACING=value][,/POLAR][,EXTRA]
KEYWORDS:
U = the X component of the two-dimensional field.
u must be a two-dimensional array.
V = the Y component of the two-dimensional field.
v must be a two-dimensional array.
X = the x-coordinate
Y = the y-coordinate
SPACING = the spacing between one arrow and the next
Higher values will produce a more "rarefied" plot
POLAR = set this keyword if the input data are in polar coordinates
(i.e. U = vr, V = vphi, x = r, y = phi) and the vector plot
should be output in cartesian coordinates.
This is useful when combined with DISPLAY, e.g.,
IDL> DISPLAY,x1=x1,x2=x2,/polar,rho
IDL> VECFIELD,v1,v2,x1,x2,/POLAR,/OVERPLOT
Note: for spherical input data (u,v)=(vr,vth) simply call
VECFIELD with x2 -> 0.5*!PI-x2 and v -> -v.
EXTRA = all the keyword accepted by velovect, such as
/overplot, length, etc...
LAST MODIFIED: Aug 27, 2012
(See vecfield.pro)
NAME: WRITE_VTK
AUTHOR: Andrea Mignone (mignone@ph.unito.it)
PURPOSE: Create a VTK data file out of a PLUTO binary file assuming
the grid information previously obtained with pload.
A scalar VTK file is created if only one argument is given;
otherwise a VTK vector file is created if 3 arguments are given.
SYNTAX: WRITE_VTK, QX[, QY, QZ][,FILENAME=string][,VARNAME=string]
[/APPEND]
ARGUMENTS:
QX = a 3D single-precision (float) binary data to be written
in the VTK file.
QY,QZ = the y- and z- components of the vector fields to be
written.
FILENAME = the name of the output VTK file. Default is "data.vtk".
VARNAME = the name of the variable appearing in the VTK header.
Default is "IDL_var".
KEYWORDS:
/APPEND: open file and append data
LAST MODIFIED:
Sept 25, 2012 by A. Mignone (mignone@ph.unito.it)
(See write_vtk.pro)