Basic+examples

=Basic librat / start operation= Librat is the library of function calls around which you can write your own code to do things such as read in and parse an object file, read in and parse camera, illumination files, waveband files and so on. However, //start// is a wrapper code around these commands which gives you access to all the basic operations, and so is the de facto tool for doing simulations. The key things required to carry out a simulation are:


 * A camera file
 * An illumination file
 * A waveband file
 * An object file - this is always assumed to be the last file on the //start// command line

Anything specific you want to do in any of these parts of the process is specified in these files. There are a limited number of additional command line options which either allow you to override a few key things in these files (the waveband file for example), or more usually are external to these things. Each of these can be passed through via the -RAT keyword. Examples are the ray tree depth (-RATm), verbose level (-RATv), waveband file (-RATsensor_wavebands) etc.

Example 1: simple object from lidar test
To run a simple lidar test, we consider a simple scene object [|test0.obj]

mtllib plants.matlib usemtl white v 0 0 0 v 0 0 1 plane -1 -2 !{ usemtl white  !{  v 0 0 1000  ell -1 30000 30000 1000  !} !}

This defines a material library [|plants.matlib] that specifies the reflectance and transmittance properties of the scene materials.

In this case, there is only a single material of type (standard reflectance material - Lambertian reflectance (and/or transmittance). The material name is  and the (ASCII) file giving the spectral reflectance function is [|refl/white.dat]. The file [|refl/white.dat] contains 2 columns: column 1 is wavelength, column 2 is reflectance for that wavelength (wavelength units are arbitrary, but we usually use nm). In this case, the file specifies:

which is a reflectance of 1.0 for any wavelength (less than or equal to an arbitrary upper limit 10000). If the file specifies transmittance as well, this is given as a third column. The lines:

specify that the material library [|plants.matlib] is to be loaded (only one library can be loaded in current versions), and that material is to be used for subsequent objects. The fields {{!{ }} and {{!} }} specifiy that a bounding box should be placed around objects contained within the brackets. This allows for efficient intersection tests in the ray tracing.

The fields denote a vertex (vector) (as in the standard wavefront format). This requires 3 numbers to be given after the giving the {x,y,z} coordinates of the vector. Note that fields can specify a location or direction vector. The fields and  specify scene objects. We will look at a fuller range of such objects later, but these two allow for a simple scene specification. is an infinite planar object. It is defined by an intersection point (location vector) and a direction vector. These vectors need to be defined before a call is made to the object, so in this case, we define as  and  as, i.e. an x-y plane at z=0. Thus means 'define a plane with  given by the previous (-1) specified vector that goes through  given by the second to last specified vector.'

is an ellipsoid object. Its description requires definition of: (i) the **base** (N.B. not the centre) of the ellipsoid ( here, meaning the previously-defined vector - in this case); (ii) the semi-axis lengths in x,y,z directions ( here). so:

is in fact a spheroid of x-y semi-axis length 30000 units (arbitrary linear units) and z-semi-axis length 1000 units: a prolate spheroid that extends from -30000 to 30000 in the x- and y-directions and from 1000 to 3000 in the z-direction.

<span style="font-size: 1.3em; margin: 0px; padding-bottom: 0px; padding-left: 0px; padding-right: 0px; padding-top: 5px;">**Example 2: running start with another simple object**
Another simple object file, [|test.obj]:

!{ mtllib plants.matlib v 0.000000 0.000000 0.000000 v 0.000000 0.000000 1.000000 usemtl soil plane -1 -2 !{ g object 0 usemtl stem v 0 0 0 v 0 0 1 cyl -1 -2 0.1 sph -1 0.2 v -1 0 1 cyl -1 -2 0.1 !} !{ clone 0 0 0 0 object 0 clone 0 1 0 90 object 0 clone -1 0 0 -90 object 0 !} !}
 * 1) define

Note the braced bounding boxes. The **mtllib** specifies the name of the material library (relative to the obj file in this case). This file points to the **locations** of the spectral information for all the materials in the object file (specified by the **usemtl** lines - soil and stem here), not the actual spectral information. The plants.matlib file might look like:

srm soil refl/price_soil.dat srm stem refl/leaf.dat

and so on. Where **srm** line (surface reflectance material) specifies a material and then the name of a file which contains the actual spectral information. Those files look like:

450.00 0.048959 460.00 0.048783 470.00 0.048640 480.00 0.048541 490.00 0.048722 500.00 0.049784 510.00 0.056469 .........

where each line is a waveband and a reflectance value. If three columns are specified the third column is assumed to be transmittance (for leaves for example). If transmittance is specified then you can (if required) define a value for t / (r + t) in the waveband file i.e. this is the ratio of transmittance to single-scattering albedo (total scattering i.e. r + t) and is assumed to be 0.5 by default. This value is used to weight the probability of transmittance to reflectance at each interaction (with this material). So if for example r >> t for a specific wavelength then a weighting of 0.5 will be inefficiently sampling as many t rays as r rays, when the r rays are far more significant. In this case, specifying t/(r+t) will weight the sampling properly.

<span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;">The object has a soil plane (material soil), defined by a "vertex" **v** which is the location

<span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;">v 0.000000 0.000000 0.000000

<span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;">and a normal

<span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;">v 0.000000 0.000000 1.000000

<span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;">i.e. normal is pointing vertically up. Then there is a #define group (g) called "object 0" is defined, which uses material stem, and is composed of a cylinder of radius 0.1, with base at 0 0 0 and tip at 0 0 1. The line

<span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;">cyl -1 -2 0.1

<span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;">specifies a cylinder with the bottom and top -1 and -2 vertex lines above this point in the file respectively. Then we have a sphere, located at the last vertex (v 0 0 1) and of radius 0.2, then another cylinder. So "cyl", "sph" and so on are the primitives that the ray tracer works with - it can also use ellipsoids (ell), closed cylinders (ccyl) and more basically, facets (f) defined by vertices e.g.:

<span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> v 0 0 0 <span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> v 0 1 0 <span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> v 1.41 1.41 0 <span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> f -3 -2 -1

<span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> defines an equilateral triangle in the x,y plane of side length 1, with 1 vertex at 0, 0, 0.

<span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;">Returning to the file above, another cylinder is specified which has its base at 1 0 1 and its tip in the same location as <span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;">the tip of the first one (the share this vertex line). This means the object is asymmetric so we can see the effect of rotating it during cloning. So although we have defined an object (**object 0**), we haven't actually invoked it yet. The lines:

<span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> !{ <span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> clone 0 0 0 0 object 0 <span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> clone 0 1 0 90 object 0 <span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> clone -1 0 0 -90 object 0 <span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> !}

<span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> define actual instances of the object "object 0". So we now clone the object 3 times, with the numbers in each case being

<span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;">x, y, z, azimuth_angle

<span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;">Running start
<span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;">You can run //start// direct in interactive mode just by typing //start// on the cmd line. Thi waits for input and if you enter an unspecified option number it will list the available options:

<span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;">pallas% start <span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;">23 <span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;">options: <span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> -2 : print PID <span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> -1 : print memory use <span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> 0 : quit <span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> 1 n s1x s1y s1z ... : set sun vectors <span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> 2 : print sun vectors <span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> 3 n b1 w1 ...i bn wn: set wavebands <span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> 4 : print wavebands <span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> 5 file.obj : read object file <span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> 6 fx fy fz dx dy dz : trace ray from f in direction d <span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> 7 : get and print materials <span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> 8 : print object information <span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> 9 : print info on materials used <span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> 10 : get and set verbosity level (0-1) <span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> 11 : get and print object bbox information <span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> 12 : <span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> 13 : same as 14 assuming filenames camera.dat light.dat <span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> 14 camera.dat light.dat : ray tracing using defined camera & illumination <span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> 15 : dont go there <span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;"> 16 cx cy cz sx sy nrows ncols rpp name : produce a height map in name

The simplest operation is to generate a s<span style="border-collapse: collapse; font-family: arial,sans-serif; line-height: normal;">imple height map (option 16)

set obj = test.obj pallas% set obj = test.obj pallas% echo 16 0 0 1000 3 3 256 256 1 $obj.zbuf | nice +19 start -v $obj

so here the option to start is 16, and we look from 0, 0, 1000 units at an area of 3 x 3 units, with an op image of 256 x 256 pixels, 1 ray per pixel and with an output image name of $obj.zbuf. The code outputs the progress of the simulation to the stderr channel. The op of the simulation is a HIPS format image (**XXXX**), a single frame of floating point values giving the height (in scene units) of the first object intercepted by a single ray in each pixel. If we look at the stats of the image using a HIPS utility called hipstats (included with the librat distribution in $BPMS/bin/x86_64/hipstats) we see:

pallas% hipstats test.obj.zbuf Frame 0 : n=65536 mean=0.097 SD=0.313 Min=0.000 Max=1.200

i.e. the max height in the scene is 1.2 units as we would hope. We can display this image using either xv (which needs /home/mdisney/bin/x86_64 to be in your $path variable), by converting it to a byte format HIPS image using another HIPS utility (also included with the librat distribution), linear. This performs a linear contrast stretch, converts the output to byte format (0-255) and then writes the result to stdout so it can be piped to xv:

pallas% linear < test.obj.zbuf | xv -

which should look like:



We can do thesame via a wrapper script (xvc), which is in /home/plewis/bpms/bin/csh so you also need this in $path i.e.

pallas% xvc test.obj.zbuf

This is a quick way to test whether the obj file works as we expect it to and looks sensible.

Converting a HIPS file to something more generic
We can use hips2pbm, a HIPS utility (also included with the librat distribution) to convert to a pbm file (portable bit map):

pallas% linear < test.obj.zbuf | hips2pbm > test.obj.zbuf.pbm

and then display the result using any number of things including ImageMagick:

pallas% display test.obj.zbuf.pbm

and we can use ImageMagick to any other format we like, via convert:

pallas% convert test.obj.zbuf.pbm test.obj.zbuf.gif

Running start using camera and light files
More often, we want to use start to simulate some type of sensor / observation characteristics. In this case we specify a camera file, a light file and the various spectral and simulation characteristics accompanying them. Then we can use option 14 to start, which expects the next two files to be the camera and light files respectively.

A camera file looks like [|camera.test.dat] (note all lines starting # are ignored - these are either comments or provide additional options which we will explore later):

camera { camera.name = "simple camera"; geometry.perspective = TRUE;
 * 1) camera type - others are spherical, planar, albedo
 * 1) default is orthographic

geometry.azimuth = 0.0; geometry.zenith = 0.0; geometry.twist = 0.0; geometry.lookAt = 0, 0, 0; geometry.boomLength = 10.; geometry.idealArea = 5; geometry.boomLength = 10.; samplingCharacteristics.nPixels = 262144; samplingCharacteristics.rpp = 4;
 * 1) geometry
 * 1) geometry.fieldOfView = FOV;

result.integral.mode = "scattering order" result.integral.format = "ascii" result.image = "op.hips" result.integral = "op.results"
 * 1) basic o/p
 * 1) result.integral.mode = "waveband"


 * 1) sampling options
 * 2) result.samplingImage = "samplingImage.hips"
 * 3) samplingPattern.OPImage = "testOPImage.hips"
 * 4) samplingPattern.sampleCountImage = "testsampleCountImage.hips"
 * 5) samplingPattern.gridMap = "testgridMap.hips"
 * 6) samplingPattern.size = 512,512;
 * 7) samplingPattern.form = "circular"
 * 8) samplingPattern.form = "gaussian"
 * 9) samplingPattern.sd = 50,50;
 * 10) samplingPattern.centre = 0,0;

}
 * 1) lidar options
 * 2) result.integral = "results_lidar.dat"
 * 3) result.integral.mode = "distance"
 * 4) lidar.nBins = 150
 * 5) lidar.binStart = 17.6
 * 6) lidar.binStep = 0.1

A light file looks like [|light.test.dat]:

camera { camera.name = "simple illumination"; geometry.azimuth = 0.0; geometry.zenith = 0.0; }

and we can run a simulation of this type as:

pallas% echo 14 camera.test.dat light.test.dat | nice +19 start -v -RATm 5 -RATsensor_wavebands wb.test.dat $obj

where -RATm 5 specifies 5 orders of scattering, the verbose level is set and the waveband file [|wb.test.dat] looks like this:

0 450 1 550 2 650 3 850

The result is a HIPS image op.hips which has 4 frames, so we can generate a fcc image as follows:

pallas% reseq 3 2 1 < op.hips |linear | hips2pbm | xv -

The actual simulated data are in the file [|op.results.direct] and look like:

1 4.373408e-02 7.572531e-02 9.658440e-02 1.624538e-01 2 2.082191e-05 6.769709e-05 1.074853e-04 4.715291e-04 3 3.809157e-06 2.112662e-05 4.297797e-05 2.893788e-04 4 2.480746e-08 2.563754e-07 6.507600e-07 1.090333e-05 5 1.032238e-09 1.880159e-08 6.054241e-08 1.903158e-06
 * 1) integral results, ordered by: scattering order (rows 5 ) waveband (cols 4 ) distance (frames 1)
 * 2) Distance (start of bin): -1.000000e+00
 * 3) Wavelength:
 * 4) 4.500000e+02 5.500000e+02 6.500000e+02 8.500000e+02

so the first 5 orders of scattering are as above (4 wb, left to right). To calculate the total scattering in each waveband we simply add up the individual interactions i.e.

pallas% grep -v # op.results.direct | gawk '{for(i=2;i<=NF;i++) r[i-1]+=$i;}END{for(i=1;i<=NF-1;i++) print r[i];}' 0.0437587

0.0758144 0.0967356 0.163228

And if we want to plot reflectance as a function of waveband we can do:

pallas% grep -v # op.results.direct | gawk '{for(i=2;i<=NF;i++) r[i-1]+=$i;}END{for(i=1;i<=NF-1;i++) print r[i];}' | paste wb.test.dat - | gawk '{print $2, $3;}' > op.results.direct.plot pallas% generate_graph op.results.direct.plot

Now let's look at a much more complex example of a single birch tree, using the file [|Birch_new1.tmp.obj]. We can edit the camera file to give a different o/p by changing the lines:

result.image = "birch.op.hips" result.integral = "birch.op.results"

we can check t

pallas% set obj = Birch_new1.tmp.obj pallas% echo 14 camera.test.dat light.test.dat | nice +19 start -v -RATm 5 -RATsensor_wavebands wb.test.dat $obj pallas% reseq 3 2 1 < birch.op.hips | linear | hips2pbm | xv -

which looks a funny colour as linear stretches the range 0 - 1 to 0 - 255 by default. If we look at the image stats we see there's not much above 0.5 even in the NIR, so we can set a range for the contrast stretch which will make better use of the colours e.g.

pallas% reseq 3 2 1 < birch.op.hips | linear 0 0.6 | hips2pbm | xv -



Building scenes: combining multiple object files with #include
We can make combine multiple object files efficiently using #include e.g. in [|test.include.obj]:

mtllib plants.matlib !{ v 0.000000 0.000000 0.000000 v 0.000000 0.000000 1.000000 usemtl soil plane -1 -2 !} clone 0 0 0 0 Birch_new1.tmp.obj 0 clone -2 0 0 55.35 Birch_new1.tmp.obj 0 !} !}
 * 1) include Birch_new1.include.obj

Note the object file Birch_new1.include.obj is slightly modified from above - we have to remove the plants.matlib line and the soil plane as we already have one of those in our global object. However now we can duplicate this plant wherever we want in the scene at little extra computational cost.



Clones can be rotated individually around the x, y, z axes by using Rx, Ry and Rz labels, followed by an angle e.g.
 * Rotations and transformations**

clone 0 0.555556 0 Ry 137.508

Clones can also make use of the povray type translation and scaling transforms (see http://graphics.stanford.edu/~cek/rayshade/raynotes/6-transf/transfor.htm for more details). In this case, the transformation is defined as:

clone 0 0 0 Transform -0.001937 0.04222 -0.09063 -0.09989 -0.004583 0 0.004154 -0.09054 -0.04226 object

which transforms (translates, rotates and scale) the specified object. The first three numbers specify the standard x, y, z clone translation (move/shift) operation. The nine numbers following the Transform label (x1 - z3) describe a 3 x 3 transformation matrix, which can provide scaling and shear.

This allows mapping of an image on to the surface of an object. To use a material map, each vertex in the [|object file] must have a local (2D) coordinate associated with it e.g.
 * Material maps**

<span style="font-family: Courier New,Courier;">mtllib [|plants.matlib] <span style="font-family: Courier New,Courier;">v -20 -20 1 <span style="font-family: Courier New,Courier;">!n -1 0 0 1 <span style="font-family: Courier New,Courier;">!local -1 0 0 <span style="font-family: Courier New,Courier;">v -20 20 1 <span style="font-family: Courier New,Courier;">!n -1 0 0 1 <span style="font-family: Courier New,Courier;">!local -1 0 1 <span style="font-family: Courier New,Courier;">v 20 20 1 <span style="font-family: Courier New,Courier;">!n -1 0 0 1 <span style="font-family: Courier New,Courier;">!local -1 1 1 <span style="font-family: Courier New,Courier;">v 20 -20 1 <span style="font-family: Courier New,Courier;">!n -1 0 1 0 <span style="font-family: Courier New,Courier;">!local -1 1 0 <span style="font-family: Courier New,Courier;">usemtl leaf <span style="font-family: Courier New,Courier;">f -2 -3 -4 <span style="font-family: Courier New,Courier;">f -1 -2 -4

Here, 4 vertices are defined, along with associated normal vector (<span style="font-family: Courier New,Courier;">!n vertexNumber dx dy dz ) and local coordinate (within the range 0 to 1). The object represented above consists of two [-2 -3 -4|triangles], forming a square, defined over the local coordinate system (0,0) to (1,1). The material named <span style="font-family: Courier New,Courier;">[|leaf] is used for these triangles. The file [|plants.matlib]contains a definition of the material <span style="font-family: Courier New,Courier;">leaf :

<span style="font-family: Courier New,Courier;">matmap leaf leaf.matmap

which means the material named leaf is a material map (matmap), which is defined in the file [|leaf.matmap].

<span style="font-family: Courier New,Courier;">material_map: plewis.hips <span style="font-family: Courier New,Courier;">attributes: <span style="font-family: Courier New,Courier;">0 material0 <span style="font-family: Courier New,Courier;">1 material1 <span style="font-family: Courier New,Courier;">... etc

In this file, the actual map is given as the file [|plewis.hips], a byte HIPL-format image: The attributes section of the file [|leaf.matmap]. associates a material name with each DN in the image [|plewis.hips]. Thus, DN 0 is associated with material0, DN 1 with material1. The material names material0, material1 etc, are defined in [|plants.matlib].