5. Examples

Now that we have a reasonable idea how NEMO is structured and used, we should be ready to go through some real examples Some of the examples below are short versions of scripts available online in one of the directories (check $NEMO/scripts and $NEMOBIN). The manual pages programs(8NEMO) and intro(1NEMO) are useful to find (and cross-reference) programs if you’re a bit lost. Each program manual should also have some references to closely related programs, in the “SEE ALSO” section.

5.1. N-body experiments

In this section we will describe how to set up an N-body experiment, run, display and analyze it. In the first example, we shall set up a head-on collision between two spherical “galaxies” and do some simple analysis.

5.1.1. Setting it up

In Filestructure we already used mkplummer to create a Plummer model; so here we shall use the program mkommod (“MaKe an Osipkov-Merritt MODel”) to make two random N-body realizations of a King model with dimensionless central potential \(W_c = 7\) and 100 particles each. The small number of particles is solely for the purpose of getting results within a reasonable time. Adjust it to whatever you can afford on your CPU and test your patience and integrator.

1% mkommod in=$NEMODAT/k7isot.dat out=tmp1 nbody=100 seed=123

These models are produced in so-called RMS-units in which the gravitational constant \(G=1\), the total mass \(M=1\), and binding energy \(E=-1/2\). In case you would like virial units (also known as N-body units) (see also: Heggie & Mathieu, \(E=-1/4\), in: The use of supercomputers in stellar dynamics ed. Hut & McMillan Springer 1987, pp.233) the models have to be rescaled using snapscale:

2% snapscale in=tmp1 out=tmp1s rscale=2 "vscale=1/sqrt(2.0)"

Note the use of the quotes in the expression, to prevent the shell to give special meaning to the parenthesis, which are shell meta characters.

The second galaxy is made in a similar way, with a different seed (the default seed=0 takes the time of the day):

3% mkommod in=$NEMODAT/k7isot.dat out=tmp2 nbody=100 seed=987

This second galaxy needs to be rescaled too, if you want virial units:

4% snapscale in=tmp2 out=tmp2s rscale=2 "vscale=1/sqrt(2.0)"

We then set up the collision by stacking the two snapshots, albeit with a relative displacement in phase space. The program snapstack was exactly written for this purpose:

5% snapstack in1=tmp1s in2=tmp2s out=i001.dat deltar=4,0,0 deltav=-1,0,0

The galaxies are initially separated by 4 unit length and approaching each other with a velocity consistent with infall from infinity (parabolic encounter). The particles assembled in the data file i001.dat are now ready to be integrated.

To look at the initials conditions we could use:

6% snapplot i001.dat xrange=-5:5 yrange=-5:5

which is displayed in Figure X

5.1.2. Integration using hackcode1

We then run the collision for 20 time units, with the standard N-body integrator based on the Barnes “hierarchical tree” algorithm:

7% hackcode1 in=i001.dat out=r001.dat tstop=20 freqout=2 \
   freq=40 eps=0.05 tol=0.7 options=mass,phase,phi > r001.log

The integration frequency relates to the integration timestep as freq=40, the softening length eps=0.05, and opening angle or tolerance tol=theta. A major output of masses, positions and potentials of all particles is done every 1/freqout = 0.5 time units, which corresponds to about 1/5 of a crossing time. The standard output of the calculation is diverted to a file r001.log for convenience. This is an (ASCII) listing, containing useful statistics of the run, such as the efficiency of the force calculation, conserved quantities etc. Some of this information is also stored in diagnostic sets in the structured binary output file r001.dat.

As an exercise, compare the output of the following two commands:

8% more r001.log
9% tsf r001.dat | more

5.1.3. Display and Initial Analysis

As described in the previous subsection, hackcode1 writes various diagnostics in the output file. A summary of conservation of energy and center-of-mass motion can be graphically displayed using snapdiagplot:

10% snapdiagplot in=r001.dat

The program snapplot displays the evolution of the particle distribution, in projection (in the yt package this is called a phase plot):

11% snapplot in=r001.dat

Depending on the actual graphics (yapp) interface snapplot was compiled with, you may have to hit the RETURN key, push a MOUSE BUTTON or just WAIT to advance from one to the next frame.

The snapplot program has a very powerful tool built into it which makes it possible to display any projection the user wants.

As an example consider:

12% snapplot in=r001.dat xvar=r yvar="x*vy-y*vx" xrange=0:10 yrange=-2:2 \

plots the angular momentum of the particles along the \(z\) axis, \(J_z = x.v_y - y.v_x\) , against their radius, \(r\), but only for the even numbered particles, (i%2==0) within a distance of 0.2 of the X-Y plane (\(-0.2<z \& z<0.2\)). Again note that some of the expressions are within quotes, to prevent the shell of giving them a special meaning.

The xvar, yvar and visib expressions are fed to the C compiler (during runtime!) and the resulting object file is then dynamically loaded into the program for execution. The expressions must contain legal C expressions and depending on their nature must return a value in the context of the program. E.g. xvar and yvar must return a real value, whereas visib must return a boolean (false/true or 0/non-0) value. This should be explained in the manual page of the corresponding programs.

In the context of snapshots, the expression can contain basic body variables which are understood to the bodytrans(3NEMO) routine. The real variables x, y, z, vx, vy, vz are the cartesian phase-space coordinates, t the time, m the mass, phi the potential, ax,ay,az the cartesian acceleration and aux some auxiliary information. The integer variables are i, the index of the particle in the snapshot (0 being the first one in the usual C tradition) and key, another spare slot.

For convenience a number of expressions have already been pre-compiled (see also Table ref{t:bodytrans}), e.g. the radius r= \(\sqrt{x^2+y^2+z^2}\) = sqrt(x*x+y*y+z*z), and velocity v = \(\sqrt{v_x^2+v_y^2+v_z^2}\) = sqrt(vx*vx+vy*vy+vz*vz). Note that r and v themselves cannot be used in expressions, only the basic body variables listed above can be used in an expression.

When you need a complex expression that has be used over and over again, it is handy to be able to store these expression under an alias for later retrieval. With the program bodytrans it is possible to save such compiled expressions object files under a new name.

Some precompiled bodytrans expressions















key (see also real version below)












(x*ax+y*ay+z*az)/sqrt(x*x+y*y+z*z) or: (rvec$cdot$avec)/$|$rvec$|$ \


















phi+0.5*(vx*vx+vy*vy+vz*vz) or: $phi$ + vvec$^2$/2 \






sqrt(sqr(x*vy-y*vx)+sqr(y*vz-z*vy)+sqr(z*vx-x*vz)) or: $|$rvec$times$vvec$|$ \



key (see also int} version above)










































\(l\) , atan2(y,x)*180/PI [-180,180]



\(b\) , atan2(z,sqrt(x*x+y*y))*180/PI [-90,90]



(-vx\sin{l} + vx\cos{l})/r



\((-vx\cos{l}\sin{b} - vy\sin{l}\sin{b}+vz\cos{b})/r\)



Aitoff projection x [-2,2] T.B.A.



Aitoff projection y [-1,1] T.B.A.



Right Ascension, converted from radians to degrees (-x)



Declination, converted from radians to degrees (y)

Here is an example of its use:

13% bodytrans expr="x*vy-y*vz" type=real file=jz

this saves the expression for the angular momentum in a real valued bodytrans expression file, btr_jz.so which can in future programs be referenced as expr=jz (whenever a real-valued bodytrans expression is required), e.g.

14% snapplot i001.dat xvar=r yvar=jz xrange=0:5

Alternatively, one can handcode a bodytrans function, compile it, and reference it locally. This is useful when you have truly complicated expressions that do not easily write themselves down on the commandline. The (x,y) AITOFF projection are an example of this. For example, consider the following code in a (local working directory) file btr_r2.c:

#include <bodytrans.h>

real btr_r2(b,t,i)
Body *b;
real t;
int  i;
    return sqrt(x*x + y*y);

By compiling this:

15% bake btr_r2.so

a shared object file btr_r2.so is created in the local directory, which could be used in any real-valued bodytrans expression:

16% snapplot i001.dat xvar=r2 yvar=jz xrange=0:5

For this your environment variable BTRPATH must have been set to include the local working directory, designated by a dot. Normally your NEMO system manager will have set the search path such that the local working directory is searched before the system one (in $NEMOOBJ/bodytrans).

5.1.4. Advanced Analysis

more to come

5.1.5. Generating models

more to come

5.1.6. Using Unix pipes

In most cases a NEMO file can be used in a pipe (usually via in=- and out=-), therefore limiting the need to write files. Here is an example of plotting the measured and expected surface brightness of a homogeneous sphere of 1,000,000 particles with unit mass and unit radius:

1 % mkconfig - 1000000 ball seed=0 |\
2     snapgrid - - nx=800 ny=800  |\
3     ccdellint - 0:1.1:0.01 inc=0 out=- |\
4     ccdprint - x=  newline=t label=x |\
5     tabmath - - '1.5/pi*sqrt(1-%1**2)' |\
6     tabplot -  1 2,3 color=2,3 line=0,0,1,1 point=2,0.1,0,0 \
7       xlab="Radius" ylab="Surface Brightness" headline="mkconfig shape=ball" \
8       yapp=ball.png/png

A few comments on the highlighted lines: In line 3 the out= keyword is not the second keyword, hence the explicit way it was written with the out=-. In line 5 the expected surface brightness expression is added as the 3rd column to the table in the pipe, then passed on for a quick and dirty plot (shown below).


5.1.7. Handling large datasets

One of NEMOs weaknesses is arguably also its strong point: programs must generally be able to fit all their data in (virtual) memory. Although programs usually free memory associated with data that is not needed anymore, there is a very clear maximum to the number of particles it can handle in a snapshot. By default a particle takes up about 100 bytes, which limits the size of a snapshots.

It may happen that your data was generated on a machine which had a lot more memory then the machine you want to analyze your data on. As long as you have the diskspace, and as long as you don’t need programs that cannot operate on data in serial mode, there is a solution to this problem. Instead of keeping all particles in one snapshot, they are stored in several snapshots of (equal number of) bodies, and as long as all snapshots have the same time and are stored back to back, most programs that can operate serially, will do this properly and know about it. Of course it’s best to split the snapshots on the machine with more memory

% snapsplit in=run1.out out=run1s.out nbody=10000

If it is just one particular program (e.g. snapgrid that needs a lot of extra memory, the following may work:

% snapsplit in=run1.out out=- nbody=1000 times=3.5 |\
    snapgrid in=- out=run1.ccd nx=1000 ny=1000 stack=t

Using tcppipe(1NEMO) one can also pipe data from the machine with large memory (A) to your workstation with smaller memory (B), as can be demonstrate with the following code snippet:

A% mkplummer - 1000 | tcppipe
B% tcppipe A | tsf -

5.2. Tables

Tables are obiquitous in astronomy in general and in NEMO specifically. Often they are simple ASCII tables, each row on a separate line and columns separated by a space or a comma (CSV). NEMO handles both, though generally output tables are space separated. NEMO’s table.5 man page goes in some more detail.

Intelligent tables (e.g. IPAC tables, the astropy ECSV format) add some header type information to annotate the table with column information (name, type, units) and sometimes even FITS-like provenance (e.g. IPAC tables). NEMO has limited ways to process these, currently loosing much of this header and provenance.

There is no good consistency in the programs if an out= keyword is needed, or if the table is written to stdout, or if out=- is the last parameter. Otherwise table programs are well suited for working the {tt pipe and filter} paradigm, as most NEMO programs are. Some examples are given below, with the caveat that the output pipe construction can be inconsistent.

Columns are typically given through the xcol= and ycol= keywords, and are 1-based. Column 0 has a special meaning, it is the row number (again 1 being the first row).

5.2.1. Manipulation

Here are some programs that manipulate tables:

  • tabgen was designed to produce tables from scratch, with some type of random values. Mostly to test performance and scalability of tables. Here a small table of 4 (rows) by 3 (columns) is written to *stdout:

% tabgen - nr=4 nc=3
0.552584 0.126911 0.520753
0.0930232 0.563683 0.258931
0.0369577 0.965763 0.634585
0.981766 0.334841 0.988963
  • tabcols selects columns from a table

  • tabrows selects rows from a table

  • tabtab takes two or more tables, and will paste (increasing the columns) or catenate (increasing the rows)

  • tabtranspose transposes the columns and rows of a table

% tabgen - nr=4 nc=3 | tabtranspose  - -
0.552584 0.0930232 0.0369577 0.981766
0.126911 0.563683 0.965763 0.334841
0.520753 0.258931 0.634585 0.988963
  • tabcomment comments out things that look like comments. Sometimes this filter is needed before another table program can be used, as they can get confused with non-numeric information. For example picking out two columns from a dirty table:

% tabcomment mytable.tab | tabcols - 2,3
  • tabcsv converts a table to some XSV, where X can be choosen from a set of characters

  • tabmath is a column calculator, but columns are referenced with a % sign, so %2 refers to column 2.. By default it will add the new columns to the old columns, unless some selection, or all, is/are removed. It can also select rows based on the values in a row, for example, only select rows where column 2 is positive.

    Here is an example of adding two columns

% tabgen -  nr=4 nc=3 | tabmath - - %1+%2 all
  • txtpar extracts parameters from a text file. The often complex operations involving Unix tools such as grep/awk/sed/cut/head/tail can often be simplified with txtpar. Here is an example :

% cat example.txt
Worst fractional energy loss dE/E = (E_t-E_0)/E_0 = 0.00146761 at T = 1.28125
QAC_STATS: - 0.039966 0.0274195 0.00185505 0.135854  0.799319 1  20

% txtpar in=example.txt expr="log(abs(%1)),log(abs(%3/%2))" format=%.3f \
         p0=Worst,1,9  p1=QAC,1,3 p2=QAC,1,4
-2.833 -0.164

5.2.2. Generating

Some other NEMO programs that produce tables:

  • snapprint tabulates choosen body variables in ascii. The output format can also be choosen.

% mkplummer - 3 seed=123 | snapprint - comment=t format=%8.5f
#       x       y       z       vx      vy      vz
-2.19118 -0.01225  0.18687 -0.21951 -0.14248  0.28165
 3.22756 -0.27674 -0.88792  0.27564  0.18469  0.11529
-1.03638  0.28899  0.70105 -0.05613 -0.04221 -0.39694
  • orblist tabulates an orbit (should be renamed to orbprint)

  • ccdprint tabulates values of an image

5.2.3. Plotting

  • tabhist creates a histogram

  • tabplot creates a line or scatterplot of one or more columns

  • tabview , tabzoom dynamic queries (cf. glueviz, topcat)

5.2.4. Analysis

  • tabtrend computes differences between successive rows

  • tabint integrate a sorted table

  • tabsmooth performs a smoothing operation on one or more columns

  • tablsqfit and tabnllsqfit do linear and non-linear fitting of functions

  • tabstat gives various statistics on selected columns

% tabgen - nr=4 nc=3 | tabstat - 1:3
npt:      4         4         4
min:      0.0369577 0.126911  0.258931
max:      0.981766  0.965763  0.988963
sum:      1.66433   1.9912    2.40323
mean:     0.416083  0.4978    0.600808
disp:     0.382993  0.311225  0.262247
skew:     0.424316  0.39325   0.250172
kurt:    -1.43956  -1.19861  -1.07587
min/sig: -0.989902 -1.19171  -1.30365
max/sig:  1.47701   1.50362   1.48011
median:   0.322804  0.449262  0.577669

5.3. Orbits

In this section we will describe how to integrate individual stellar orbits, display and analyze them. Be aware that although 3D orbits can be computed the number of utilities to analyze them is rather limited.

Orbits are normally stored in datafile (see also orbit(5NEMO)), and a close conceptual relationship exists between a (single-particle type) snapshot and an orbit: an orbit is an ordered series of phase-space coordinates whereas a snapshot is a series of particles with no particular order, but all at the same time.

Since orbits will be computed in an analytical potential, we assume for the remainder of this section that you have familiarized yourself with how to supply potentials to orbit integrator programs. They all share the same triple potname=, potpars=, potfile= keyword interface, as described in Section ref{s:potential}. Many examples of the tricky potpars= keyword are given in Appendix ref{a:potential}.

5.3.1. Initializing

There are a few programs with which orbits can be initialized:

  • mkorbit is the most straightforward program. You can give simply give it all 6 phase space coordinates, and an orbit file consisting of this one point is generated. It is also possible to give the potential in which the particle is to move, and 5 phase space coordinates plus the energy, or even 4 phase space coordinates and the energy plus the total angular momentum or angular momentum along the Z axis (for axisymmetric systems).

Let’s start with an example of creating a simple orbit by itself with no associated potential.

% mkorbit out=orb1 x=1 y=0 z=0 vx=0 vy=0.2 vz=0
### Warning [mkorbit]: Potential potname= not used; set etot=0.0
pos: 1.000000 0.000000 0.000000
vel: 0.000000 0.200000 0.000000
etot: 0.000000

% tsf orb1
char History[59] "mkorbit out=orb1 x=1 y=0 z=0 vx=0 vy=0.2 vz=0 VERSION=3.2b"
set Orbit
  set Parameters
    int Ndim 03
    double Mass 1.00000
    double IOM[3] 0.00000 0.200000 0.00000
    int Nsteps 01
  set Potential
  set Path
    double TimePath[1] 0.00000
    double PhasePath[1][2][3] 1.00000 0.00000 0.00000 0.00000 0.200000 0.00000
  • perorb is a program that for given initial conditions (similar to the ones described in {tt mkorbit} above) attempts to calculate periodic orbits in that potential. The output file will be a file with one (or more) orbits. This is a bit of an advanced program, and will not be covered here.

  • stoo is a program that can take a particle position from a snapshot, and turn it into an orbit. For example, sampling some initial conditions from the positions of stars in a Plummer sphere, we could use the following small (bash) shell code to find some statistical properties of this selected set of orbitsfootnote{For the careful reader: mkplummer and potname=plummer actually have different units, and as such this experiment is not properly set up.}

mkplummer out=p100 nbody=p100
for i in $(nemoinp 0:100:10); do
    stoo in=p100 out=orb$i ibody=$i
    orbint orb$i orb$i.out 10000 0.01 10000 potname=plummer
    orbstat orb$i.out

The reverse program, otos turns an orbit into a snapshot, and may come in handy since the snapshot package has far more advanced analysis programs.

5.3.2. Integration

  • orbint integrates orbits from given initial conditions. If the input orbit has more than 1 step, the last step is taken as the initial conditions. Although the potname=, potpars=, potfile= keywords can be given, if the input orbit contains…

  • perorb finds periodic orbits, and stores a full period which should close the orbit. This program finds periodic orbits in the XY plane (i.e. currently it will only find 2D orbits) by searching for the centers of invariant curves in the surface of section.

  • henyey also finds periodic orbits, but uses Henyey’s methodfootnote{see also van Albada & Sanders, (1982, MNRAS, 201, 303)}. This program has however not been released to the public version of NEMO, and in fact it seems the source code was lost.

5.3.3. Display

  • orbplot is the only orbit plotting program we currently have. For more sophisticated display tabplot and/or {tt snapplot} would have to be used after transforming the data. Also {tt snapplot} uses the powerful bodytrans expression parser to plot arbitrary body related expressions, although orbplot can handle both x, y, z and vx, vy, vz for the xvar= and yvar= keywords. An example of the output of orbplot is given in Figure ref{f:orbit1}.

5.3.4. Analysis

  • orbstat is an example of a simple program that reads orbits, and displays statistics of it’s 2D (x-y-) coordinates: maximum extent, as well as statistics of the angylar momentum. This program is not suited for 3D orbits yet.

% orbint orb1 orb1.long
% orbstat orb1.out
# T     E       x_max   y_max   u_max   v_max   j_mean  j_sigma
1000 -0.687107 1 0.999958 0.746764 0.746611 0.2 3.83111e-09
  • orbfour performs a variety of fourier analysis on the coordinates

% orbint orb1 orb1.long 100000 0.01 10000 10 plummer
INIPOTENTIAL Plummer: [3d version]
Pattern speed=0
0.000000 0.020000 -0.707107     -0.6871067811865
100.000000 0.277794 -0.964901     -0.6871067811856
200.010000 0.020912 -0.708019     -0.6871067812165
300.020000 0.271222 -0.958329     -0.6871067812194
400.030000 0.023376 -0.710483     -0.6871067812465
500.040000 0.259253 -0.946360     -0.6871067812551
600.050000 0.027415 -0.714522     -0.6871067812765
700.060000 0.242979 -0.930086     -0.6871067812904
800.070000 0.033056 -0.720163     -0.6871067813065
900.080000 0.223694 -0.910801     -0.6871067813241
Energy conservation: 2.00138e-10
% orbfour orb1.long amode=t
<R> N A0 A1 A2 A3 A4 B1 B2 B3 B4
1 10001 0.000360461 0.334714     0.000150399 -0.000472581 -0.000158864
                   -0.000667155  0.000228086 -0.000725406  0.000103029

% orbfour orb1.long amode=f
<R> N C0 C1 P1 C2 P2 C3 P3 C4 P4
1 10001 0.000360461
        0.334715      -0.114202
        0.000273209   56.5992
        0.000865763 -123.083
        0.000189349  147.035
  • orbsos computes surface of section coordinates. Since this program does not plot, but produces a simple ascii table, you can pipe the output into tabplot:

% orbsos orb1.long y | tabplot - 3 4  xlab=Y ylab=VY
% orbsos orb1.long x | tabplot - 3 4  xlab=X ylab=VX

will plot either a Y-VY or X-VX surface of section.

  • orbdim computes the dimensionality of an orbit, i.e. how many integrals of motions it has. Although it requires very long integration times to accurately compute this, it is completely automatic, and does not require an analysis like that for a surface of section (which is also graphic). It is based on an interesting paper by Carnevali & Santangelo (1984, ApJ 281 473-476).

  • otos transforms an orbit back into a snapshot, thereby giving you the much richer set of analysis tools that are available for snapshot.

5.4. Potential(tmp)

Here we list some of the standard potentials available in NEMO, in a variety of units, so not always \(G=1\)!

Recall that most NEMO program use the keywords potname= for the identifying name, potpars= for an optional list of parameters and potfile= for an optional text string,for example for potentials that need some kind of text file. The parameters listed in potpars= will always have as first parameter the pattern speed in cases where rotating potentials are used. A Plummer potential with mass 10 and core radius 5 would be hence be supplied as: potname=plummer potpars=0,10,5. The plummer potential ignored the potfile= keyword.

plummer Plummer potential (BT, pp.42, eq. 2.47)

\[\Phi = - { M \over { {(r_c^2 + r^2)}^{1/2} } }\]
  • \(\Omega_p\) : Pattern Speed (always the first parameter in potpars=)

  • \(M\) : Total mass (clearly G=1 here)

  • \(r_c\) : Core radius

potname=plummer potpars= \(\Omega_p\), \(M\), \(r_c\)

5.5. Images

5.5.1. Initializing Images

There are a few programs with which images can be initialized:

  • ccdgen is the simplest program to generate images from scratch, much like tabgen for tables. Here is an example creating and visualizing an image

% ccdgen bar.ccd  object=bar spar=1,10,0.5,30 size=128,128
% ccdplot bar.ccd
% ccdfits bar.ccd  bar.fits
% ds9 bar.fits
  • ccdmath is the most straightforward program. Here is an example of creating an image from scratch:

% ccdmath out=ccd1 fie=%x+%y size=2,4
Generating a map from scratch

% tsf ccd1
set Image
  set Parameters
    int Nx 2
    int Ny 4
    int Nz 1
    double Xmin 0.00000
    double Ymin 0.00000
    double Zmin 0.00000
    double Dx 1.00000
    double Dy 1.00000
    double Dz 1.00000
    double Xrefpix 0.00000
    double Yrefpix 0.00000
    double Zrefpix 0.00000
    double MapMin -4.00000
    double MapMax 0.00000
    int BeamType 0
    double Beamx 0.00000
    double Beamy 0.00000
    double Beamz 0.00000
    double Restfreq 0.00000
    double VLSR 0.00000
    double Time 0.00000
    char Storage[5] "CDef"
    int Axis 1
  set Map
    double MapValues[2][4] 0.00000 1.00000 2.00000 3.00000 1.00000
       2.00000 3.00000 4.00000

% ccdprint ccd1 x= y= label=x,y
 Y\X 0 1

3  3 4
2  2 3
1  1 2
0  0 1
  • snapgrid converts a snapshot to an image.

  • fitsccd converts a FITS file to an image. The inverse of this, ccdfits also exists. Storage model of FITS vs. NEMO’s image format is currently not optimal, and thus for very large images a CPU penalty is incurred.

nx,ny ->    data[nx][ny]

e.g.  ccdmath out=ccd1   nx=10 ny=5
      gives       double MapValues[10][5]

ccdmath "" - %x 3,2 | tsf - margin=100 | grep MapVal
      MapValues[3][2] -2.00000 -2.00000 -1.00000 -1.00000 0.00000 0.00000
ccdmath "" - %y 3,2 | tsf - margin=100 | grep MapVal
      MapValues[3][2] -2.00000 -1.00000 -2.00000 -1.00000 -2.00000 -1.00000

5.5.2. Galactic Velocity Fields

As an example, a special section is devoted here to the analysis of galactic velocity fields.footnote{In this example shell variables such as r=$(nemoinp 0:60) have been replaced with the more portable macro files like @tmp.r. Although the example uses 0:60 and works fine in the shell the example was used under, increasing the number to 256 would fail because of overflowing the maximum characters allowed on the commandline}

The following programs are available:

    ccdvel          create a model velocity field, from scratch
    rotcur          tilted ring model velocity field fitting
    rotcurshape     annulus rotation curve shape fitting to a velocity field
    ccdmath         perform math on images, or use math to create images
    ccdplot         plot (contour/greyscale) an image
    ccdprint        print out pixel values in an imamge

% nemoinp 0:60 > tmp.r
% tabmath tmp.r - "100*%1/(20+%1)" all > tmp.v
% ccdvel out=map1.vel rad=@tmp.r vrot=@tmp.v pa=30 inc=60
% rotcurshape in=map1.vel radii=0,60 pa=30 inc=60 vsys=0 units=arcsec,1 \
              rotcur1=core1,100,20,1,1 tab=-

% ccdmath out=map0.vel fie=0 size=128,128
% rotcurshape map0.vel 0,40 30 45 0 blank=-999 resid=map2.vel \
              rotcur1=plummer,200,10,0,0 fixed=all units=arcsec,1

Since rotcurshape computes a residual velocity field, one can easily create nice model velocity fields from any selected shape by fitting a rotation curve shape to a velocity field of all 0s and keeping all parameters fixed to the requested values:

% ccdmath out=map0.vel fie=0 size=128,128
% rotcurshape map0.vel 0,40 30 45 0 blank=-999 resid=map.vel \
           rotcur1=plummer,200,10,0,0 fixed=all units=arcsec,1
% ccdplot map.vel -100:100:10 blankval=0 cmode=1

5.6. Exchanging data


examples/Exchanging Data