source: Sophya/trunk/SophyaLib/Manual/sophya.tex@ 1594

Last change on this file since 1594 was 1436, checked in by ansari, 25 years ago

derniere (?) MAJ doc sophya.tex avant tag - Reza 9/3/2001

File size: 49.6 KB
Line 
1\documentclass[twoside,11pt]{article}
2% Package standard : Utilisation de caracteres accentues, mode francais et graphique
3\usepackage{url}
4\usepackage[latin1]{inputenc}
5\usepackage[T1]{fontenc}
6\usepackage[english]{babel}
7\usepackage{graphicx}
8% package a mettre pour faire du pdf
9\usepackage{palatino}
10
11% Extension de symboles mathematiques
12\usepackage{amssymb}
13
14% Definition pour Docs Sophya
15\usepackage{defsophya}
16
17% Constitution d'index
18\usepackage{makeidx}
19\makeindex
20
21\begin{document}
22
23\begin{titlepage}
24% The title page - top of the page with the title of the paper
25\titrehp{Sophya \\ An overview }
26% Authors list
27\auteurs{
28R. Ansari & ansari@lal.in2p3.fr \\
29E. Aubourg & aubourg@hep.saclay.cea.fr \\
30G. Le Meur & lemeur@lal.in2p3.fr \\
31C. Magneville & cmv@hep.saclay.cea.fr \\
32S. Henrot-Versille & versille@in2p3.fr
33}
34% \auteursall
35% The title page - bottom of the page with the paper number
36\vspace{1cm}
37\begin{center}
38{\bf \Large Sophya Version: 1.1 (V\_Fev2001) \\
39Document revision 1.0 }
40\end{center}
41\titrebp{1}
42\end{titlepage}
43
44\tableofcontents
45
46\newpage
47
48\section{Introduction}
49
50{\bf SOPHYA} ({\bf SO}ftware for {\bf PHY}sics {\bf A}nalysis)
51is a collection of C++ classes designed for numerical and
52physics analysis software development. Our goal is to provide
53easy to use, yet powerful classes which can be used by scientists.
54We have decided to use as much as possible available
55numerical analysis libraries, encapsulating them whenever
56possible. Although some the modules in SOPHYA have been
57designed with the specific goal of providing the general framework for
58the Planck-HFI data processing software, most of the
59packages presented here have a more general scope than the
60CMB analysis and Planck mission problem.
61\par
62\vspace*{2mm}
63This documents
64presents only a brief overview of the class library,
65mainly from the user's point of view. A more complete description
66can be found in the reference manual, available from our
67web site: {\bf http://www.sophya.org}.
68\par
69\vspace*{2mm}
70The source directory tree
71\footnote{ CVS: cvsserver.lal.in2p3.fr:/exp/eros/CVSPlanck}
72is organised into a number of modules.
73
74\begin{itemize}
75\item[] {\bf Mgr/} Scripts for code management,
76makefile generation and software installation
77\item[] {\bf SysTools/} General architecture support classes such
78as {\tt PPersist, NDataBlock<T>}, and few utility classes
79({\tt DataCard, DVList} \ldots). Presently, this module contains
80also classes implementing interfaces to OS specific services, such
81as shared object and dynamic link handling. This module will be
82separated in three modules, one for the general architecture
83support, one for the utility classes, and one for the OS specific
84services.
85\item[] {\bf TArray/} template numerical arrays, vectors and matrices \\
86({\tt PixelMap<T> SphericalMap<T>} \ldots)
87\item[] {\bf NTools/} Some standard numerical analysis tools
88(linear, and non linear parameter fitting, FFT, \ldots)
89\item[] {\bf HiStats/} Histogram-ming and data set handling classes \\
90({\tt Histo Histo2D NTuple XNTuple} \ldots)
91\end{itemize}
92
93The modules listed below are more tightly related to the
94CMB (Cosmic Microwave Background) data analysis problem:
95\begin{itemize}
96\item[] {\bf SkyMap/} Local and full sky maps, and few geometry
97handling utility classes. \\
98({\tt PixelMap<T>, LocalMap<T>, SphericalMap<T>, \ldots})
99\item[] {\bf SkyT/}
100classes for spectral emission and detector frequency response modelling \\
101({\tt SpectralResponse, RadSpectra, BlackBody} \ldots)
102\item[] {\bf Samba/} Spherical harmonic analysis, noise generators \ldots
103\end{itemize}
104
105The following modules contain the interface classes with
106external libraries:
107\begin{itemize}
108\item[] {\bf FitsIOServer/} Classes for handling file input-output
109in FITS format using the cfitsio library.
110\item[] {\bf LinAlg/} Interface with Lapack linear algebra package
111\item[] {\bf IFFTW/} Interface with FFTW package (libfftw.a)
112\end{itemize}
113
114Other modules:
115\begin{itemize}
116\item[] {\bf Tests/} Simple test programs
117\item[] {\bf PrgUtil/} Various utility programs (runcxx, scanppf, scanfits, \ldots)
118\item[] {\bf PMixer/} skymixer and related programs
119\item[] {\bf ProgPI/} interactive analysis tool - It should be noted that
120this module uses the SOPHYA class library and is based on {\bf PI}
121which is a C++ library defining a complete GUI program
122architecture. An additional module (PIext) define the interactive
123analysis program framework and the interfaces with the objects
124in SOPHYA. The {\bf PI/} \footnote{the PI package documentation
125is available from {\bf http://www.lal.in2p3.fr/recherche/eros/PeidaDoc/} }
126and {\bf PIext/} modules are not currently part
127of the SOPHYA CVS structure.
128\end{itemize}
129
130Modules containing examples and demo programs:
131\begin{itemize}
132\item[] {\bf Examples/} Sample SOPHYA codes and example programs and
133makefiles (auto\_makefile and ex\_makefile).
134\item[] {\bf DemoPIApp/} Sample exripts and programs for (s)piapp
135interactive analysis tools.
136\end{itemize}
137\newpage
138
139\section{Using Sophya}
140Basic usage of Sophya classes are described in in the following sections.
141Complete Sophya documentation can be found at our web site: \\
142{\bf http://www.sophya.org}.
143
144\subsection{Environment variables}
145Two environment variables {\bf DPCBASEREP} and {\bf EROSCXX} are used
146to define the path where the Sophya libraries and executable are installed.
147{\bf DPCBASEREP} defines the base directory path and {\bf EROSCXX} the
148name of the C++ compiler. The complete path is built using {\bf DPCBASEREP},
149the operating system name (as obtained by the {\tt uname} command), and
150the compiler name. In the example below, we show the complete path
151for a {\tt Linux} system, using the GNU g++ compiler:
152
153\begin{itemize}
154\item \$DPCBASEREP/Include : Include (.h) files
155\item \$DPCBASEREP/Linux-g++/Libs : Path for the archive libraries (.a)
156\item \$DPCBASEREP/Linux-g++/ShLibs : Shared library path (.so)
157\item \$DPCBASEREP/Linux-g++/Exec : Executable file path
158\end{itemize}
159
160In order to use the shared libraries, the {\bf LD\_LIBRARY\_PATH} variable
161should contain the Sophya shared library path
162({\tt \$DPCBASEREP/Linux-g++/ShLibs } when using g++ compiler on Linux)
163
164For modules using external libraries, the {\bf EXTLIBDIR}
165environment variable should contain the path to these libraries
166and corresponding include files.
167C-FitsIO anf FFTW include files should be accessible through: \\
168{\tt \$EXTLIBDIR/Include/FitsIO } \\
169{\tt \$EXTLIBDIR/Include/FFTW } \\
170The corresponding libraries are expected to be found in: \\
171{\tt \$EXTLIBDIR/Linux-g++/Libs} \\
172
173\subsection{User makefiles}
174The file {\tt \$DPCBASEREP/Include/MakefileUser.h} defines the compilation
175flags and the list of Sophya libraries. It should be included in the
176user's makefile. The default compilation rules assumes that the object (.o)
177and executable files would be put in the following directories: \\
178{\tt \$HOME/`uname`-\$EROSCXX/Objs} \\
179{\tt \$HOME/`uname`-\$EROSCXX/Exec}.
180In the case of a {\tt Linux} system and using {\tt g++} as the C++ compiler,
181these two directories would be translated to \\
182{\tt \$HOME/Linux-g++/Objs} and {\tt \$HOME/Linux-g++/Exec}.
183The GNU make program should be used.
184\par
185The file {\tt Examples/auto\_makefile} defines the rules to compile
186a given source program, and link it against the Sophya libraries to produce
187an executable. The example below shows the steps to compile a program named
188{\tt trivial.cc }.
189\begin{verbatim}
190csh> cp Examples/auto_makefile makefile
191csh> cp Examples/ex1.cc trivial.cc
192csh> make trivial
193\end{verbatim}
194This command should compile the {\tt trivial.cc} file,
195and link it against the sophya libraries.
196The file {\tt Examples/ex\_makefile} provides another
197example makefile.
198
199\subsection{the runcxx program}
200\index{runcxx}
201{\bf runcxx} is a simple program which can be used to compile, link
202and run simple C++ programs. It handles the creation of a
203complete program file, containing the basic set C++ include files,
204the necessary include files for SOPHYA SysTools, TArray, HiStats
205and NTools modules, and the main program with exception handling.
206Other Sophya modules can be included using the {\tt -import} flag.
207Use of additional include files can be specified using the
208{\tt -inc} flag.
209\begin{verbatim}
210csh> runcxx -h
211SOPHYA Version 1.1 Revision 0 (V_Fev2001) -- Feb 28 2001 11:19:17 cxx
212 runcxx : compiling and running of a piece of C++ code
213 Usage: runcxx [-compopt CompileOptions] [-linkopt LinkOptions]
214 [-tmpdir TmpDirectory] [-f C++CodeFileName]
215 [-inc includefile] [-inc includefile ...]
216 [-import modulename] [-import modulename ...]
217 [-uarg UserArg1 UserArg2 ...]
218 if no file name is specified, read from standard input
219 modulenames: SkyMap, Samba, SkyT, FitsIOServer, LinAlg, IFFTW
220\end{verbatim}
221Most examples in this manual can be tested using runcxx. The
222example below shows how to compile, link and run a sample
223code.
224\begin{verbatim}
225// File example.icc
226Matrix a(3,3);
227a = IdentityMatrix(1.);
228cout << a ;
229// Executing this sample code
230csh> runcxx -f example.icc
231\end{verbatim}
232
233\newpage
234
235\section{Copy constructor and assignment operator}
236In C++, objects can be copied by assignment or by initialisation.
237Copying by initialisation corresponds to creating an object and
238initialising its value through the copy constructor.
239The copy constructor has its first argument as a reference, or
240const reference to the object's class type. It can have
241more arguments, if default values are provided.
242Copying by assignment applies to an existing object and
243is performed through the assignment operator (=).
244The copy constructor implements this for identical type objects:
245\begin{verbatim}
246class MyObject {
247public:
248 MyObject(); // Default constructor
249 MyObject(MyObject const & a); // Copy constructor
250 MyObject & operator = (MyObject const & a) // Assignment operator
251}
252\end{verbatim}
253The copy constructors play an important role, as they are
254called when class objects are passed by value,
255returned by value, or thrown as an exception.
256\begin{verbatim}
257// A function declaration with an argument of type MyObject,
258// passed by value, and returning a MyObject
259MyObject f(MyObject x)
260{
261 MyObject r;
262 ...
263 return(r); // Copy constructor is called here
264}
265// Calling the function :
266MyObject a;
267f(a); // Copy constructor called for a
268\end{verbatim}
269It should be noted that the C++ syntax is ambiguous for the
270assignment operator. {\tt MyObject x; x=y; } and
271{\tt MyObject x=y;} have different meaning.
272\begin{verbatim}
273MyObject a; // default constructor call
274MyObject b(a); // copy constructor call
275MyObject bb = a; // identical to bb(a) : copy constructor call
276MyObject c; // default constructor call
277c = a; // assignment operator call
278\end{verbatim}
279
280As a general rule in SOPHYA, objects which implements
281reference sharing on their data members have a copy constructor
282which shares the data, while the assignment operator copies or
283duplicate the data.
284
285\newpage
286\section{Module SysTools}
287
288{\bf SysTools} contains utility classes such as {\tt DataCards} or
289{\tt DVlist}, an hierarchy of exception classes for Sophya, a template
290class {\tcls{NDataBlock}} for handling reference counting on numerical
291arrays, as well as classes providing the services for implementing simple
292serialisation.
293\vspace*{5mm}
294
295\subsection{SOPHYA persistence}
296\index{PPersist} \index{PInPersist} \index{POutPersist}
297\begin{figure}[hbt]
298\dclsa{PPersist}
299\dclsbb{PIOPersist}{PInPersist}
300\dclsb{POutPersist}
301\caption{partial class diagram for classes handling persistence in Sophya}
302\end{figure}
303A simple persistence mechanism is defined in SOPHYA. Its main
304features are:
305\begin{itemize}
306\item[] Portable file format, containing the description of the data structures
307and object hierarchy. \\
308{\bf PPF} {\bf P}ortable {\bf P}ersistence file {\bf F}ormat.
309\item[] Handling of read/write for multiply referenced objects.
310\item[] All write operations are carried using sequential access only. This
311holds also for read operations, unless positional tags are used.
312SOPHYA persistence services can thus be used to transfer objects
313through network links.
314\item[] The serialisation (reading/writing) for objects for a given class
315is implemented through a handler object. The handler class inherits
316from {\tt PPersist} class.
317\item[] A run time registration mechanism is used in conjunction with
318RTTI (Run Time Type Identification) for identifying handler classes
319when reading {\bf PInPersist} streams, or for associating handlers
320with data objects {\bf AnyDataObject} for write operations.
321\end{itemize}
322A complete description of SOPHYA persistence mechanism and guidelines
323for writing delegate classes for handling object persistence is beyond
324the scope of this document. The most useful methods for using Sophya
325persistence are listed below:
326\begin{itemize}
327\item[] {\tt POutPersist::PutObject(AnyDataObj \& o)} \\
328Writes the data object {\bf o} to the output stream.
329\item[] {\tt POutPersist::PutObject(AnyDataObj \& o, string tagname)} \\
330Writes the data object {\bf o} to the output stream, associated with an
331identification tag {\bf tagname}.
332\item[] {\tt PInPersist::GetObject(AnyDataObj \& o)} \\
333Reads the next object in stream into {\bf o}. An exception is
334generated for incompatible object types.
335\item[] {\tt PInPersist::GetObject(AnyDataObj \& o, string tagname)} \\
336Reads the object associated with the tag {\bf tagname} into {\bf o}.
337An exception is generated for incompatible object types.
338\end{itemize}
339The operators {\tt operator << (POutPersist ...) } and
340{\tt operator >> (PInPersist ...) } are often overloaded
341to perform {\tt PutObject()} and {\tt GetObject()} operations,
342as illustrated in the example below:
343\begin{verbatim}
344// Creating and filling a histogram
345Histo hw(0.,10.,100);
346...
347// Writing histogram to a PPF stream
348POutPersist os("hw.ppf");
349os << hw;
350// Reading a histogram from a PPF stream
351PInPersist is("hr.ppf");
352is >> hr;
353\end{verbatim}
354
355The {\bf scanppf} program can be used to list the content of a PPF file.
356\index{scanppf}
357\begin{verbatim}
358csh> scanppf -h
359SOPHYA Version 1.1 Revision 0 (V_Fev2001) -- Feb 28 2001 11:19:17 cxx
360 Usage: scanppf filename [s/n/a0/a1/a2/a3]
361 s[=default} : Sequential reading of objects
362 n : Object reading at NameTags
363 a0...a3 : Tag List with PInPersist.AnalyseTags(0...3)
364\end{verbatim}
365
366
367\subsection{\tcls{NDataBlock}}
368\index{\tcls{NDataBlock}}
369\begin{figure}[hbt]
370\dclsbb{AnyDataObj}{\tcls{NDataBlock}}
371\dclsbb{PPersist}{\tcls{FIO\_NDataBlock}}
372\end{figure}
373The {\bf \tcls{NDataBlock}} is designed to handle reference counting
374and sharing of memory blocs (contiguous arrays) for numerical data
375types. Initialisation, resizing, basic arithmetic operations, as
376well as persistence handling services are provided.
377The persistence handler class ({\tt \tcls{FIO\_NDataBlock}}) insures
378that a single copy of data is written for multiply referenced objects,
379and the data is shared among objects when reading.
380\par
381The example below shows writing of NDataBlock objects through the
382use of overloaded operator $ << $ :
383\begin{verbatim}
384#include "fiondblock.h"
385// ...
386POutPersist pos("aa.ppf");
387NDataBlock<r_4> rdb(40);
388rdb = 567.89;
389pos << rdb;
390// We can also use the PutObject method
391NDataBlock<int_4> idb(20);
392idb = 123;
393pos.PutObject(idb);
394\end{verbatim}
395The following sample programs show the reading of the created PPF file :
396\begin{verbatim}
397PInPersist pis("aa.ppf");
398NDataBlock<r_4> rdb;
399pis >> rdb;
400cout << rdb;
401NDataBlock<int_4> idb;
402cout << idb;
403\end{verbatim}
404
405\subsection{Using DVList}
406\index{DVList} \index{MuTyV}
407\begin{figure}[hbt]
408\dclsbb{AnyDataObj}{DVList}
409\dclsbb{PPersist}{\tclsc{ObjFileIO}{DVList}}
410\end{figure}
411The {\bf DVList} class objects can be used to create and manage list
412of values, associated with names. A list of pairs of (MuTyV, name(string))
413is maintained by DVList objects. {\bf MuTyV} is a simple class
414capable of holding string, integer, float or complex values,
415providing easy conversion methods between these objects.
416\begin{verbatim}
417// Using MuTyV objects
418MuTyV s("hello"); // string type value
419MuTyV x;
420x = "3.14159626"; // string type value, ASCII representation for Pi
421double d = x; // x converted to double = 3.141596
422x = 314; // x contains the integer value = 314
423// Using DVList
424DVList dvl;
425dvl("Pi") = 3.14159626; // float value, named Pi
426dvl("Log2") = 0.30102999; // float value, named Log2
427dvl("FileName") = "myfile.fits"; // string value, named myfile.fits
428// Printing DVList object
429cout << dvl;
430\end{verbatim}
431
432\subsection{Using DataCards}
433\index{DataCards}
434The {\bf DataCards} class can be used to read parameters from a file.
435Each line in the file starting with \@ defines a set of values
436associated with a keyword. In the example below, we read the
437parameters corresponding with the keyword {\tt SIZE} from the
438file {\tt ex.d}. We suppose that {\tt ex.d} contains the line: \\
439{\tt @SIZE 400 250} \\
440\begin{verbatim}
441#include "datacards.h"
442// ...
443// Initialising DataCards object dc from file ex.d
444DataCards dc( "ex.d" );
445// Getting the first and second parameters for keyword size
446// We define a default value 100
447int size_x = dc.IParam("SIZE", 0, 100);
448int size_y = dc.IParam("SIZE", 1, 100);
449cout << " size_x= " << size_x << " size_y= " << size_y << endl;
450\end{verbatim}
451
452\subsection{Dynamic linker}
453\index{PDynLinkMgr}
454The class {\bf PDynLinkMgr} can be used for managing shared libraries
455at run time. The example below shows the run time linking of a function:\\
456{\tt extern "C" { void myfunc(); } } \\
457\begin{verbatim}
458#include "pdlmgr.h"
459// ...
460string soname = "mylib.so";
461string funcname = "myfunc";
462PDynLinkMgr dyl(soname);
463DlFunction f = dyl.GetFunction(funcname);
464if (f != NULL) {
465// Calling the function
466 f();
467}
468\end{verbatim}
469
470\subsection{CxxCompilerLinker class}
471\index{CxxCompilerLinker}
472This class provides the services to compile C++ code and building
473shared libraries, using the same compiler and options which have
474been used to create the SOPHYA shared library.
475The sample program below illustrates using this class to build
476the shared library (myfunc.so) from the source file myfunc.cc :
477\begin{verbatim}
478#include "cxxcmplnk.h"
479// ...
480string flnm = "myfunc.cc";
481string oname, soname;
482int rc;
483CxxCompilerLinker cxx;
484// The Compile method provides a default object file name
485rc = cxx.Compile(flnm, oname);
486if (rc != 0 ) { // Error when compiling ... }
487// The BuildSO method provides a default shared object file name
488rc = cxx.BuildSO(oname, soname);
489if (rc != 0 ) { // Error when creating shared object ... }
490\end{verbatim}
491
492\newpage
493\section{Module TArray}
494\index{\tcls{TArray}}
495{\bf TArray} module contains template classes for handling standard
496operations on numerical arrays. Using the class {\tt \tcls{TArray} },
497it is possible to create and manipulate up to 5-dimension numerical
498arrays {\tt (int, float, double, complex, \ldots)}. The include
499file {\tt array.h} declares all the classes and definitions
500in module TArray. {\bf Array} is a typedef for arrays
501with double precision floating value elements. \\
502{\tt typedef TArray$<$r\_8$>$ Array ; }
503
504\begin{figure}[hbt]
505\dclsccc{AnyDataObj}{BaseArray}{\tcls{TArray}}
506\dclsbb{PPersist}{\tcls{FIO\_TArray}}
507\end{figure}
508
509
510\subsection{Using arrays}
511\index{Sequence} \index{RandomSequence} \index{RegularSequence}
512\index{EnumeratedSequence}
513The example below shows basic usage of arrays, creation, initialisation
514and arithmetic operations. Different kind of {\bf Sequence} objects
515can be used for initialising arrays.
516
517\begin{figure}[hbt]
518\dclsbb{Sequence}{RandomSequence}
519\dclsb{RegularSequence}
520\dclsb{EnumeratedSequence}
521\end{figure}
522
523The example below shows basic usage of arrays:
524\index{\tcls{TArray}}
525\begin{verbatim}
526// Creating and initialising a 1-D array of integers
527TArray<int> ia(5);
528EnumeratedSequence es;
529es = 24, 35, 46, 57, 68;
530ia = es;
531cout << "Array<int> ia = " << ia;
532// 2-D array of floats
533TArray<r_4> b(6,4), c(6,4);
534// Initializing b with a constant
535b = 2.71828;
536// Filling c with random numbers
537c = RandomSequence();
538// Arithmetic operations
539TArray<r_4> d = b+0.3f*c;
540cout << "Array<float> d = " << d;
541\end{verbatim}
542
543The copy constructor shares the array data, while the assignment operator
544copies the array elements, as illustrated in the following example:
545\begin{verbatim}
546TArray<int> a1(4,3);
547a1 = RegularSequence(0,2);
548// Array a2 and a1 shares their data
549TArray<int> a2(a1);
550// a3 and a1 have the same size and identical elements
551TArray<int> a3;
552a3 = a1;
553// Changing one of the a2 elements
554a2(1,1,0) = 555;
555// a1(1,1) is also changed to 555, but not a3(1,1)
556cout << "Array<int> a1 = " << a1;
557cout << "Array<int> a3 = " << a3;
558\end{verbatim}
559
560\subsection{Matrices and vectors}
561\index{\tcls{TMatrix}} \index{\tcls{TVector}}
562\begin{figure}[hbt]
563\dclsccc{\tcls{TArray}}{\tcls{TMatrix}}{\tcls{TVector}}
564\end{figure}
565Vectors and matrices are 2 dimensional arrays. The array size
566along one dimension is equal 1 for vectors. Column vectors
567have {\tt NCols() = 1} and row vectors have {\tt NRows() = 1}.
568Mathematical expressions involving matrices and vectors can easily
569be translated into C++ code using {\tt TMatrix} and
570{\tt TVector} objects. {\bf Matrix} and {\bf Vector} are
571typedefs for double precision float matrices and vectors.
572The operator {\bf *} beteween matrices is redefined to
573perform matrix multiplication. One can then write: \\
574\begin{verbatim}
575 // We create a row vector
576 Vector v(1000, BaseArray::RowVector);
577 // Initialize values with a random sequence
578 v = RandomSequence();
579 // Compute the vector length (norm)
580 double norm = (v*v.Transpose()).toScalar();
581 cout << "Norm(v) = " << norm << endl;
582\end{verbatim}
583
584This module contains basic array and matrix operations
585such as the Gauss matrix inversion algorithm
586which can be used to solve linear systems, as illustrated by the
587example below:
588\begin{verbatim}
589#include "sopemtx.h"
590// ...
591// Creation of a random 5x5 matrix
592Matrix A(5,5);
593A = RandomSequence(RandomSequence::Flat);
594Vector X0(5);
595X0 = RandomSequence(RandomSequence::Gaussian);
596// Computing B = A*X0
597Vector B = A*X0;
598// Solving the system A*X = B
599Vector X;
600LinSolve(A, B, X);
601// Checking the result
602Vector diff = X-X0;
603cout << "X-X0= " << diff ;
604double min,max;
605diff.MinMax(min, max);
606cout << " Min(X-X0) = " << min << " Max(X-X0) = " << max << endl;
607\end{verbatim}
608
609\subsection{Working with sub-arrays and Ranges}
610\index{Range}
611A powerful mechanism is included in array classes for working with
612sub-arrays. The class {\bf Range} can be used to specify range of array
613indexes in any of the array dimensions. Any regularly spaced index
614range can be specified, using the {\tt start} and {\tt end} index
615and an optional step (or stride). It is also possible to specify
616the {\tt start} index and the number of elements.
617In the following example, a simple low-pas filter, on a one
618dimensional stream (Vector) has been written using sub-arrays:
619
620\begin{verbatim}
621// Input Vector containing a noisy periodic signal
622 Vector in(1024), out(1024);
623 in = RandomSequence(RandomSequence::Gaussian, 0., 1.);
624 for(int kk=0; kk<in.Size(); kk++)
625 in(kk) += 2*sin(kk*0.05);
626// Compute the output vector by a simple low pass filter
627 Vector out(1024);
628 int w = 2;
629 for(int k=w; k<in.Size()-w; k++)
630 out(k) = in(Range(k-w, k+w).Sum()/(2.*w+1.);
631\end{verbatim}
632
633\subsection{Memory organisation}
634{\tt \tcls{TArray} } can handle numerical arrays with various memory
635organisation, as long as the spacing (steps) along each axis is
636regular. The five axis are labeled X,Y,Z,T,U. The examples below
637illustrates the memory location for a 2-dimensional, $N_x=4 \times N_y=3$.
638The first index is along the X axis and the second index along the Y axis.
639\begin{verbatim}
640 | (0,0) (0,1) (0,2) (0,3) |
641 | (1,0) (1,1) (1,2) (1,3) |
642 | (2,0) (2,1) (2,2) (2,3) |
643\end{verbatim}
644In the first case, the array is completely packed
645($Step_X=1, Step_Y=N_X=4$), with zero offset,
646while in the second case, $Step_X=2, Step_Y=10, Offset=10$:
647\begin{verbatim}
648 | 0 1 2 3 | | 10 12 14 16 |
649Ex1 | 4 5 6 7 | Ex2 | 20 22 24 26 |
650 | 8 9 10 11 | | 30 32 34 36 |
651\end{verbatim}
652
653For matrices and vectors, an optional argument ({\tt MemoryMapping})
654can be used to select the memory mapping, where two basic schemes
655are available: \\
656{\tt CMemoryMapping} and {\tt FortranMemoryMapping}. \\
657In the case where {\tt CMemoryMapping} is used, a given matrix line
658is packed in memory, while the columns are packed when
659{\tt FortranMemoryMapping} is used. The first index when addressing
660the matrix elements (line number index) runs along
661the Y-axis if {\tt CMemoryMapping} is used, and along the X-axis
662in the case of {\tt FortranMemoryMapping}.
663Arithmetic operations between matrices
664with different memory organisation is allowed as long as
665the two matrices have the same sizes (Number of rows and columns).
666The following code example and the corresponding output illustrates
667these two memory mappings. The {\tt \tcls{TMatrix}::TransposeSelf() }
668method changes effectively the matrix memory mapping, which is also
669the case of {\tt \tcls{TMatrix}::Transpose() } method without argument.
670
671\begin{verbatim}
672TArray<r_4> X(4,2);
673X = RegularSequence(1,1);
674cout << "Array X= " << X << endl;
675TMatrix<r_4> X_C(X, true, BaseArray::CMemoryMapping);
676cout << "Matrix X_C (CMemoryMapping) = " << X_C << endl;
677TMatrix<r_4> X_F(X, true, BaseArray::FortranMemoryMapping);
678cout << "Matrix X_F (FortranMemoryMapping) = " << X_F << endl;
679\end{verbatim}
680This code would produce the following output (X\_F = Transpose(X\_C)) :
681\begin{verbatim}
682Array X=
683--- TArray<f> ND=2 SizeX*Y*...= 4x2 ---
6841, 2, 3, 4
6855, 6, 7, 8
686
687Matrix X_C (CMemoryMapping) =
688--- TMatrix<f>(NRows=2, NCols=4) ND=2 SizeX*Y*...= 4x2 ---
6891, 2, 3, 4
6905, 6, 7, 8
691
692Matrix X_F (FortranMemoryMapping) =
693--- TMatrix<f>(NRows=4, NCols=2) ND=2 SizeX*Y*...= 4x2 ---
6941, 5
6952, 6
6963, 7
6974, 8
698\end{verbatim}
699
700\newpage
701
702\section{Module HiStats}
703\begin{figure}[hbt]
704\dclsccc{AnyDataObj}{Histo}{HProf}
705\dclsbb{AnyDataObj}{Histo2D}
706\dclsbb{AnyDataObj}{Ntuple}
707\dclsb{XNtuple}
708\caption{partial class diagram for histograms and ntuples}
709\end{figure}
710
711{\bf HiStats} contains classes for creating, filling, printing and
712doing various operations on one or two dimensional histograms
713{\tt Histo} and {\tt Histo2D} as well as profile histograms {\tt HProf}. \\
714This module also contains {\tt NTuple} and {\tt XNTuple} which are
715more or less the same that the binary FITS tables.
716
717\subsection{1D Histograms}
718\index{Histo}
719For 1D histograms, various numerical methods are provided such as
720computing means and sigmas, finding maxima, fitting, rebinning,
721integrating \dots \\
722The example below shows creating and filling a one dimensional histogram
723of 100 bins from $-5.$ to $+5.$ to create a Gaussian normal distribution
724with errors~:
725\begin{verbatim}
726#include "histos.h"
727// ...
728Histo H(-0.5,0.5,100);
729H.Errors();
730for(int i=0;i<25000;i++) {
731 double x = NorRand();
732 H.Add(x);
733}
734H.Print(80);
735\end{verbatim}
736
737\subsection{2D Histograms}
738\index{Histo2D}
739Much of these operations are also valid for 2D histograms. 1D projection
740or slices can be set~:
741\begin{verbatim}
742#include "histos2.h"
743// ...
744Histo2D H2(-1.,1.,100,0.,60.,50);
745H2.SetProjX(); // create the 1D histo for X projection
746H2.SetBandX(25.,35.); // create 1D histo projection for 25.<y<35.
747H2.SetBandX(35.,45.); // create 1D histo projection for 35.<y<45.
748H2.SetBandX(40.,55.); // create 1D histo projection for 40.<y<55.
749//... fill H2 with what ever you want
750H2.Print();
751Histo *hx = H2.HProjX();
752 hx->Print(80);
753Histo *hbx2 = HBandX(1); // Get the second X band (35.<y<45.)
754 hbx2->Print(80);
755\end{verbatim}
756
757\subsection{Profile Histograms}
758\index{HProf}
759Profiles histograms {\bf HProf} contains the mean and the
760sigma of the distribution
761of the values filled in each bin. The sigma can be changed to
762the error on the mean. When filled, the profile histogram looks
763like a 1D histogram and much of the operations that can be done on 1D histo
764may be applied onto profile histograms.
765
766\subsection{Data tables (tuples)}
767\index{NTuple}
768NTuple are memory resident tables of 32 bits floating values (float).
769They are arranged in columns. Each line is often called an event.
770These objects are frequently used to analyze data.
771Graphicals tools (spiapp) can plot a column against an other one
772with respect to various selection cuts. \\
773Here is an example of creation and filling~:
774\begin{verbatim}
775#include "ntuple.h"
776#include "srandgen.h"
777// ...
778char* nament[4] = {"i","x","y","ey"};
779r_4 xnt[4];
780NTuple NT(4,nament);
781for(i=0;i<5000;i++) {
782 xnt[0] = i+1;
783 xnt[1] = 5.*drandpm1(); // a random value between -5 and +5
784 xnt[2] = 100.*exp(-0.5*xnt[1]*xnt[1]) + 1.;
785 xnt[3] = sqrt(xnt[2]);
786 xnt[2] += xnt[3] * NorRand(); // add a random gaussian error
787 NT.Fill(xnt);
788}
789\end{verbatim}
790
791XNTuple are sophisticated NTuple : they accept various types
792of column values (double,float,int,string,...) and can handle
793very large data sets, through swap space on disk.
794\index{XNTuple}
795In the sample code below we show how to create a XNTuple
796object with four columns (double, double, int, string).
797Several entries (lines) are then appended to the table,
798which is saved to a PPF file.
799\begin{verbatim}
800#include "xntuple.h"
801// ...
802char * names[4] = {"X", "X2", "XInt","XStr"};
803// XNTuple (Table) creation with 4 columns, of integer,
804// double(2) and string type
805XNTuple xnt(2,0,1,1, names);
806// Filling the NTuple
807r_8 xd[2];
808int_4 xi[2];
809char xss[2][32];
810char * xs[2] = {xss[0], xss[1]} ;
811for(int i=0; i<50; i++) {
812 xi[0] = i; xd[0] = i+0.5; xd[1] = xd[0]*xd[0];
813 sprintf(xs[0],"X=%g", xd[0]);
814 xnt.Fill(xd, NULL, xi, xs);
815}
816// Printing table info
817cout << xnt ;
818// Saving object into a PPF file
819POutPersist po("xnt.ppf");
820po << xnt ;
821\end{verbatim}
822
823\subsection{Writing, viewing \dots }
824
825All these objects have been design to be written to or read from a persistent file.
826The following example shows how to write the previously created objects
827into such a file~:
828\begin{verbatim}
829//-- Writing
830{
831char *fileout = "myfile.ppf";
832string tag;
833POutPersist outppf(fileout);
834tag = "H"; outppf.PutObject(H,tag);
835tag = "H2"; outppf.PutObject(H2,tag);
836tag = "NT"; outppf.PutObject(NT,tag);
837} // closing ``}'' destroy ``outppf'' and automatically close the file !
838\end{verbatim}
839
840Sophya graphical tools (spiapp) can automatically display and operate
841all these objects.
842
843\newpage
844\section{Module NTools}
845
846This module provides elementary numerical tools for numerical integration,
847fitting, sorting and ODE solving. FFTs are also provided (Mayer,FFTPack).
848
849\subsection{Fitting}
850\index{Fitting} \index{Minimisation}
851Fitting is done with two classes {\tt GeneralFit} and {\tt GeneralFitData}
852and is based on the Levenberg-Marquardt method.
853\index{GeneralFit} \index{GeneralFitData}
854GeneralFitData is a class which provide a description of the data
855to be fitted. GeneralFit is the fitter class. Parametrized functions
856can be given as classes which inherit {\tt GeneralFunction}
857or as simple C functions. Classes of pre-defined functions are provided
858(see files fct1dfit.h and fct2dfit.h). The user interface is very close
859from that of the CERN {\tt Minuit} fitter.
860Number of objects (Histo, HProf \dots ) are interfaced with GeneralFit
861and can be easily fitted. \\
862Here is a very simple example for fitting the previously created NTuple
863with a Gaussian~:
864\begin{verbatim}
865#include "fct1dfit.h"
866// ...
867
868// Read from ppf file
869NTuple nt;
870{
871PInPersist pis("myfile.ppf");
872string tag = "NT"; pis.GetObject(nt,tag);
873}
874
875// Fill GeneralData
876GeneralData mGdata(nt.NEntry());
877for(int i=0; i<nt.NEntry(); i++)
878 mGdata.AddData1(xnt[1],xnt[2],xnt[3]); // Fill x, y and error on y
879mGData.PrintStatus();
880
881// Function for fitting : y = f(x) + noise
882Gauss1DPol mFunction; // gaussian + constant
883
884// Prepare for fit
885GeneralFit mFit(&mFunction); // create a fitter for the choosen function
886mFit.SetData(&mGData); // connect data to the fitter
887
888// Set and initialise the parameters (that's non-linear fitting!)
889// (num par, name, guess start, step, [limits min and max])
890mFit.SetParam(0,"high",90.,1..);
891mFit.SetParam(1,"xcenter",0.05,0.01);
892mFit.SetParam(2,"sigma",sig,0.05,0.01,10.);
893 // Give limits to avoid division by zero
894mFit.SetParam(3,"constant",0.,1.);
895
896// Fit and print result
897int rcfit = mFit.Fit();
898mFit.PrintFit();
899if(rcfit>0) {)
900 cout<<"Reduce_Chisquare = "<<mFit.GetChi2Red()
901 <<" nstep="<<mFit.GetNStep()<<" rc="<<rcfit<<endl;
902} else {
903 cout<<"Fit_Error, rc = "<<rcfit<<" nstep="<<mFit.GetNStep()<<endl;
904 mFit.PrintFitErr(rcfit);
905}
906
907// Get the result for further use
908TVector<r_8> ParResult = mFit.GetParm();
909cout<<ParResult;
910\end{verbatim}
911
912Much more usefull possibilities and detailed informations might be found
913in the HTML pages of the Sophya manual.
914
915\subsection{Polynomial}
916\index{Polynomial} \index{Poly} \index{Poly2}
917Polynomials of 1 or 2 variables are supported ({\tt Poly} and {\tt Poly2}).
918Various operations are supported~:
919\begin{itemize}
920\item elementary operations between polynomials $(+,-,*,/) $
921\item setting or getting coefficients
922\item computing the value of the polynomial for a given value
923 of the variable(s),
924\item derivating
925\item computing roots (degre 1 or 2)
926\item fitting the polynomial to vectors of data.
927\end{itemize}
928Here is an example of polynomial fitting~:
929\begin{verbatim}
930#include "poly.h"
931// ...
932Poly pol(2);
933pol[0] = 100.; pol[1] = 0.; pol[2] = 0.01; // Setting coefficients
934TVector<r_8> x(100);
935TVector<r_8> y(100);
936TVector<r_8> ey(100);
937for(int i=0;i<100;i++) {
938 x(i) = i;
939 ey(i) = 10.;
940 y(i) = pol((double) i) + ey(i)*NorRand();
941 ey(i) *= ey(i)
942}
943
944TVector<r_8> errcoef;
945Poly polfit;
946polfit.Fit(x,y,ey,2,errcoef);
947
948cout<<"Fit Result"<<polfit<<endl;
949cout<<"Errors :"<<errcoef;
950\end{verbatim}
951
952Similar operations can be done on polynomials with 2 variables.
953
954\subsection{Integration, Differential equations}
955\index{Integration}
956The NTools module provide also simple classes for numerical integration
957of functions and differential equations.
958\begin{figure}[hbt]
959\dclsbb{Integrator}{GLInteg}
960\dclsb{TrpzInteg}
961\end{figure}
962
963\index{GLInteg} \index{TrpzInteg}
964{\bf GLInteg} implements the integration through Gauss-Legendre method
965and {\bf TrpzInteg} implements trapeze integration. For {\bf TrpzInteg},
966number of steps specify the number of trapeze, and integration step,
967their width.
968The sample code below illustrates the use of TrpzInteg class:
969\begin{verbatim}
970#include "integ.h"
971// ......................................................
972// Function to be integrated
973double myf(double x)
974{
975// Simple a x + b x^2 (a=2 b=3)
976return (x*(2.+3.*x));
977}
978// ......................................................
979
980// Compute Integral(myf, 2., 5.) between xmin=2., xmax=5.
981TrpzInteg trpz(myf, 2., 5.);
982// We specify an integration step
983trpz.DX(0.01);
984// The integral can be computed as trpz.Value()
985double myf_integral = trpz.Value();
986// We could have used the cast operator :
987cout << "Integral[myf, 2., 5.]= " << (double)trpz << endl;
988// Limits can be specified through ValueBetween() method
989cout << "Integral[myf, 0., 4.]= " << trpz.ValueBetween(0.,4.) << endl;
990\end{verbatim}
991
992\subsection{Fourier transform (FFT)}
993\index{FFT} \index{FFTPackServer}
994An abstract interface for performing FFT operations is defined by the
995{\bf FFTServerInterface} class. The {\bf FFTPackSever} class implements
996one dimensional FFT, on real and complex data. FFTPackServer uses an
997adapted and extended version of FFTPack (available from netlib),
998translated in C, and can operate on single and double precision
999({\tt float, double}) data.
1000
1001The sample code below illustrates the use of FFTServers:
1002\begin{verbatim}
1003#include "fftpserver.h"
1004 // ...
1005TVector<r_8> in(32);
1006TVector< complex<r_8> > out;
1007in = RandomSequence();
1008FFTPackServer ffts;
1009ffts.setNormalize(true); // To have normalized transforms
1010cout << " FFTServer info string= " << ffts.getInfo() << endl;
1011cout << "in= " << in << endl;
1012cout << " Calling ffts.FFTForward(in, out) : " << endl;
1013ffts.FFTForward(in, out);
1014cout << "out= " << out << endl;
1015\end{verbatim}
1016
1017\newpage
1018\section{Module SkyMap}
1019\begin{figure}[hbt]
1020\dclsbb{AnyDataObj}{PixelMap}
1021\dclsccc{PixelMap}{Sphericalmap}{SphereHEALPix}
1022\dclsc{SphereThetaPhi}
1023\dclsb{LocalMap}
1024\caption{partial class diagram for pixelization classes in Sophya}
1025\end{figure}
1026The {\bf SkyMap} module provides classes for creating, filling, reading pixelized spherical and 2D-maps. The types of values stored in pixels can be int, float, double , complex etc. according to the specialization of the template type.
1027\subsection {Spherical maps}
1028There are two kinds of spherical maps according pixelization algorithms. SphereHEALPix represents spheres pixelized following the HEALPIix algorithm (E. Hivon, K. Gorski)
1029\footnote{see the HEALPix Homepage: http://www.eso.org/kgorski/healpix/ }
1030, SphereThetaPhi represents spheres pixelized following an algorithm developed at LAL-ORSAY. The example below shows creating and filling of a SphereHEALPix with nside = 8 (it will be 12*8*8= 768 pixels) :
1031\index{\tcls{SphereHEALPix}}
1032\index{\tcls{SphereThetaPhi}}
1033
1034\begin{verbatim}
1035#include "spherehealpix.h"
1036// ...
1037SphereHEALPix<double> sph(8);
1038for (int k=0; k< sph.NbPixels(); k++) sph(k) = (double)(10*k);
1039\end{verbatim}
1040
1041SphereThetaPhi is used in a similar way with an argument representing number of slices in theta (Euler angle) for an hemisphere.
1042\index{\tcls{SphereThetaPhi}}
1043
1044\subsection {Local maps}
1045\index{\tcls{LocalMap}}
1046A local map is a 2 dimensional array, with i as column index and j as row index. The map is supposed to lie on a plan tangent to the celestial sphere in a point whose coordinates are (x0,y0) on the local map and (theta0, phi0) on the sphere. The range of the map is defined by two values of angles covered respectively by all the pixels in x direction and all the pixels in y direction (SetSize()). Default value of (x0, y0) is middle of the map, center of pixel(nx/2, ny/2).
1047
1048Internally, a map is first defined within this reference plane and tranported until the point (theta0, phi0) in such a way that both axes are kept parallel to meridian and parallel lines of the sphere. The user can define its own map with axes rotated with respect to reference axes (this rotation is characterized by angle between the local parallel line and the wanted x-axis-- method SetOrigin(...))
1049
1050The example below shows creating and filling of a LocalMap with 4 columns and 5 rows. The origin is set to default. The map covers a sphere portion defined by two angles of 30. degrees (methods \textit{SetOrigin()} and \textit{SetSize()} must be called in order to completely define the map).
1051\begin{verbatim}
1052#include "localmap.h"
1053//..............
1054 LocalMap<r_4> locmap(4,5);
1055 for (int k=0; k<locmap.NbPixels();k++) locmap(k)=10.*k;
1056 locmap.SetOrigin();
1057 locmap.SetSize(30.,30.);
1058\end{verbatim}
1059
1060\subsection{Writing, viewing \dots }
1061
1062All these objects have been design to be written to or read from a persistant file.
1063The following example shows how to write the previously created objects
1064into such a file~:
1065\begin{verbatim}
1066//-- Writing
1067
1068#include "fiospherehealpix.h"
1069//................
1070
1071char *fileout = "myfile.ppf";
1072POutPersist outppf(fileout);
1073FIO_SphereHEALPix<r_8> outsph(sph);
1074outsph.Write(outppf);
1075FIO_LocalMap<r_8> outloc(locmap);
1076outloc.Write(outppf);
1077// It is also possible to use the << operator
1078POutPersist os("sph.ppf");
1079os << outsph;
1080os << outloc;
1081\end{verbatim}
1082
1083Sophya graphical tools (spiapp) can automatically display and operate
1084all these objects.
1085
1086\newpage
1087\section{Module Samba}
1088\index{Spherical Harmonics}
1089\index{SphericalTransformServer}
1090The module provides several classes for spherical harmonic analysis. The main class is \textit{SphericalTranformServer}. It contains methods for analysis and synthesis of spherical maps. The following example fills a vector of Cl's, generate a spherical map from these Cl's. This map is analysed back to Cl's...
1091\begin{verbatim}
1092#include "skymap.h"
1093#include "samba.h"
1094....................
1095
1096// Generate input spectra a + b* l + c * gaussienne(l, 50, 20)
1097int lmax = 92;
1098Vector clin(lmax);
1099for(int l=0; l<lmax; l++) {
1100 double xx = (l-50.)/10.;
1101 clin(l) = 1.e-2 -1.e-4*l + 0.1*exp(-xx*xx);
1102}
1103
1104// Compute map from spectra
1105SphericalTransformServer<r_8> ylmserver;
1106int m = 128; // HealPix pixelisation parameter
1107SphereHEALPix<r_8> map(m);
1108ylmserver.GenerateFromCl(map, m, clin, 0.);
1109// Compute power spectrum from map
1110Vector clout = ylmserver.DecomposeToCl(map, lmax, 0.);
1111\end{verbatim}
1112
1113\newpage
1114\section{Module SkyT}
1115\index{RadSpectra} \index{SpectralResponse}
1116The SkyT module is composed of two types of classes:
1117\begin{itemize}
1118\item{} one which corresponds to an emission spectrum of
1119radiation, which is called RadSpectra
1120\item{} one which corresponds to the spectral response
1121of a given detector (i.e. corresponding to a detector
1122filter in a given frequency domain), which is called
1123SpectralResponse.
1124\end{itemize}
1125\begin{figure}[hbt]
1126\dclsbb{RadSpectra}{RadSpectraVec}
1127\dclsb{BlackBody}
1128\dclsccc{AnyDataObj}{SpectralResponse}{SpecRespVec}
1129\dclsc{GaussianFilter}
1130\caption{partial class for SkyT module}
1131\end{figure}
1132
1133\begin{verbatim}
1134#include "skyt.h"
1135// ....
1136// Compute the flux from a blackbody at 2.73 K through a square filter
1137BlackBody myBB(2.73);
1138// We define a square filter from 100 - 200 GHz
1139SquareFilter mySF(100,200);
1140// Compute the filtered integrated flux :
1141double flux = myBB.filteredIntegratedFlux(mySF);
1142\end{verbatim}
1143
1144A more detailed description of SkyT module can be found in:
1145{\it The SkyMixer (SkyT and PMixer modules) - Sophya Note No 2. }
1146available also from Sophya Web site.
1147
1148\newpage
1149\section{Module FitsIOServer}
1150\begin{figure}[hbt]
1151\dclsbb{FitsFile}{FitsInFile}
1152\dclsb{FitsOutFile}
1153\end{figure}
1154\index{FITS} \index{FitsInFile} \index{FitsOutFile}
1155This module provides classes for handling file input-output in FITS format using the cfitsio library. It works like the SOPHYA persistence (see Module SysTools), using delegate objects, but its design is simpler. The following example writes a matrix (see module TArray) and a spherical map (see module SkyMap) on a FITS file and reads back from FITS file and creates new objects :
1156\begin{verbatim}
1157#include "spherehealpix.h"
1158#include "fitsspherehealpix.h"
1159#include "fitstarray.h"
1160#include "tmatrix.h"
1161//...........................
1162
1163int m=...;
1164SphereHEALPix<r_8> sph(m);
1165................
1166int dim1=...;
1167int dim2=...;
1168TMatrix<r_8> mat(dim1,dim2);
1169............
1170
1171FITS_SphereHEALPix<r_8> sph_temp(sph);
1172FITS_TArray<r_8> mat_temp(mat);
1173// writing
1174
1175FitsOutFile os("myfile.fits");
1176sph_temp.Write(os);
1177mat_temp.Write(os);
1178
1179// reading
1180FitsInFile is("myfile.fits");
1181sph_temp.Read(is);
1182mat_temp.Read(is);
1183SphereHEALPix<r_8> new_sph=(SphereHEALPix<r_8>)sph_temp;
1184TMatrix<r_8> new_mat=(TMatrix<r_8>)mat_temp;
1185................
1186
1187\end{verbatim}
1188
1189The operators {\tt operator << (FitsOutFile ...)} and
1190{\tt operator >> (FitsInFile ...)} are defined in order
1191to facilitate the FITS file operations:
1192\begin{verbatim}
1193// Writing an array object to a FITS file
1194#include "fitstarray.h"
1195FitsOutFile fio("arr.fits");
1196Matrix m(20,30);
1197m = 12345.;
1198fio << m;
1199// .....
1200// Reading a binary table to a XNTuple
1201#include "fitsxntuple.h"
1202XNTuple xn;
1203FitsInFile fii("table.fits");
1204fii >> xn;
1205\end{verbatim}
1206
1207
1208\newpage
1209\section{LinAlg and IFFTW modules}
1210An interface to use LAPACK library (available from {\tt http://www.netlib.org})
1211is implemented by the {\bf LapackServer} class, in module LinAlg.
1212\index{LapackServer}.
1213The sample code below shows how to use SVD (Singular Value Decomposition)
1214through LapackServer:
1215\begin{verbatim}
1216#include "intflapack.h"
1217// ...
1218// Use FortranMemoryMapping as default
1219BaseArray::SetDefaultMemoryMapping(BaseArray::FortranMemoryMapping);
1220// Create an fill the arrays A and its copy AA
1221int n = 20;
1222Matrix A(n , n), AA;
1223A = RandomSequence(RandomSequence::Gaussian, 0., 4.);
1224AA = A; // AA is a copy of A
1225// Compute the SVD decomposition
1226Vector S; // Vector of singular values
1227Matrix U, VT;
1228LapackServer<r_8> lpks;
1229lpks.SVD(AA, S, U, VT);
1230// We create a diagonal matrix using S
1231Matrix SM(n, n);
1232for(int k=0; k<n; k++) SM(k,k) = S(k);
1233// Check the result : A = U*SM*VT
1234Matrix diff = U*(SM*VT) - A;
1235double min, max;
1236diff.MinMax(min, max);
1237cout << " Min/Max difference Matrix (?=0) , Min= " << min
1238 << " Max= " << max << endl;
1239\end{verbatim}
1240
1241\index{FFTWServer}
1242The {\bf FFTWServer} class (in module FFTW) implements FFTServerInterface class
1243methods, for one dimensional and multi-dimensional Fourier
1244transforms on double precision data using the FFTW package
1245(available from {\tt http://www.fftw.org}).
1246
1247\newpage
1248\section{Building and installing Sophya}
1249Presently, the Sophya library has been tested with the following
1250compiler/platform pairs:
1251
1252\begin{center}
1253\begin{tabular}{ll}
1254Compaq/DEC OSF1 & cxx (6.0 , 6.2) \\
1255Compaq/DEC OSF1 & g++ (2.95) \\
1256Linux & g++ (2.91 , 2.95) \\
1257Linux & KCC (3.4) \\
1258Solaris & g++ (2.95) \\
1259SGI IRIX64 & CC (7.3) \\
1260\end{tabular}
1261\end{center}
1262
1263Some of the modules in the Sophya package uses external libraries. The
1264{\bf FitsIOServer} is the example of such a module, where the {\tt libcfitsio.a}
1265is used.
1266The build procedure expects to find the include files and the libraries in: \\
1267{\tt \$EXTLIBDIR/Include/FitsIO } \\
1268{\tt \$EXTLIBDIR/`uname`-\$EROSCXX/Libs} \\
1269
1270The object files from a given Sophya module are grouped in an archive library
1271with the module's name ({\tt libmodulename.a}). All Sophya modules
1272 are grouped in a single shared library ({\tt libsophya.so}), while the
1273modules with reference to external libraries are grouped in
1274({\tt libextsophya.so}). The {\bf PI} and {\bf PIext} modules are
1275grouped in ({\tt libPI.so}).
1276
1277The environment variables {\bf DPCDEVREP}, {\bf EXTLIBDIR} and {\bf EROSCXX}
1278must be defined in order to install the Sophya package.
1279In the example below, we assume that we want to install Sophya from a
1280released (tagged) version in the source directory {\tt \$SRC} in the
1281{\tt /usr/local/Sophya} diretory, using {\tt g++}. We assume that
1282the external libraries directory tree has been set up in
1283{\tt /usr/local/ExtLibs/}. \\[3mm]
1284\centerline{
1285 \rule{20mm}{0.5mm}
1286 {\bf \large the use of GNU make is mandatory}
1287 \rule{20mm}{0.5mm} }
1288
1289\vspace*{3mm}
1290\begin{verbatim}
1291# We select our C++ compiler
1292csh> setenv EROSCXX g++
1293# Setup the build directory
1294csh> mkdir /usr/local/Sophya/
1295csh> setenv DPCDEVREP /usr/local/Sophya/
1296csh> setenv EXTLIBDIR /usr/local/ExtLibs/
1297# Use the top level makefile in Mgr/
1298csh> cd \$SRC
1299csh> cp Mgr/Makefile Makefile
1300# Step 1: Create the directory tree and copy the include files (.h)
1301csh> make depend
1302# Step 2: Compile the modules without external library reference
1303csh> make libs
1304# Step 3: Compile the modules WITH external library reference (optional)
1305csh> make extlibs
1306# Step 4: Build libsophya.so
1307csh> make slb
1308# Step 5: Build libextsophya.so (optional)
1309csh> make slbext
1310# Step 6: Compile the PI and PIext modules (optional)
1311csh> make PI
1312# Step 7: Build the corresponding shared library libPI.so (optional)
1313csh> make slbpi
1314\end{verbatim}
1315
1316To compile all modules and build the shared libraries, it is possible
1317to use:
1318\begin{verbatim}
1319# Step 2,3,6
1320csh> make all
1321# Step 4,5,7
1322csh> make slball
1323\end{verbatim}
1324
1325At this step, all libraries should have been made. Programs using
1326Sophya libraries can now be built:
1327\begin{verbatim}
1328# To compile test programs
1329csh> cd Tests
1330csh> make arrt ...
1331csh> cd ..
1332# To compile other programs, for example from the PMixer module
1333csh> cd PMixer
1334csh> make
1335csh> cd ..
1336# To build (s)piapp (libPI.so is needed)
1337csh> cd ProgPI
1338csh> make
1339csh> cd ..
1340\end{verbatim}
1341
1342\subsection{Mgr module}
1343This module contains scripts which can be used for generating the
1344makefiles for each module.
1345\begin{itemize}
1346\item {\bf Makefile} Top level Makefile for building the libraries.
1347\item {\bf Makefile.h} contains the definition of compilation flags for the
1348different compilers and systems. This file is used for building the
1349library and generating {\bf MakefileUser.h} (to be included in makefiles).
1350\item {\bf Makefile.slb} contains the rules for building shared libraries
1351for the different compilers and systems. (to be included in makefiles)
1352\item {\bf crerep\_sophya} c-shell script for creating the directory tree
1353under {\tt \$DPCBASEREP} and {\tt \$DPCDEVREP}
1354\item {\bf install\_sophya} c-shell script for installing the Sophya package.
1355Usually from {\tt \$DPCDEVREP} to {\tt \$DPCBASEREP}
1356\item {\bf mkmflien} c-shell script for making symbolic links or copying
1357include files to {\tt \$DPCDEVREP/Include} or {\tt \$DPCBASEREP/Include}
1358\item {\bf mkmf} c-shell script for generating module makefiles and the
1359top level makefile (named GNUmakefile)
1360\item {\bf mkmflib} c-shell script for generating each library module
1361makefile (named GNUmakefile)
1362\item {\bf mkmfprog} c-shell script for generating makefile for a module
1363containing the source for executable programs (named GNUmakefile).
1364\item {\bf mkmfPI} c-shell script for generating makefile for PI and PIext
1365modules (named GNUmakefile)
1366\item {\bf libdirs} List of Sophya modules without reference to external
1367libraries.
1368\item {\bf extlibdirs} List of Sophya modules with reference to external
1369libraries.
1370
1371\end{itemize}
1372
1373\newpage
1374\appendix
1375\section{SOPHYA Exceptions}
1376\index{Exception classes} \index{PThrowable} \index{PError} \index{PException}
1377SOPHYA library defines a set of exceptions which are used
1378for signalling error conditions. The figure below shows a partial
1379class diagram for exception classes in SOPHYA.
1380\begin{figure}[hbt]
1381\dclsbb{PThrowable}{PError}
1382\dclscc{PError}{AllocationError}
1383\dclscc{PError}{NullPtrError}
1384\dclscc{PError}{ForbiddenError}
1385\dclscc{PError}{AssertionFailedError}
1386\dclsbb{PThrowable}{PException}
1387\dclscc{PException}{IOExc}
1388\dclscc{PException}{SzMismatchError}
1389\dclscc{PException}{RangeCheckError}
1390\dclscc{PException}{ParmError}
1391\dclscc{PException}{TypeMismatchExc}
1392\dclscc{PException}{MathExc}
1393\dclscc{PException}{CaughtSignalExc}
1394\caption{partial class diagram for exception handling in Sophya}
1395\end{figure}
1396
1397For simple programs, it is a good practice to handle
1398the exceptions at least at high level, in the {\tt main()} function.
1399The example below shows the exception handling and the usage
1400of Sophya persistence.
1401
1402\input{ex1.inc}
1403
1404
1405\newpage
1406\addcontentsline{toc}{section}{Index}
1407\printindex
1408\end{document}
1409
Note: See TracBrowser for help on using the repository browser.