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

Last change on this file since 3037 was 3037, checked in by ansari, 19 years ago

doc SophyaOverview complete en vue du tag V2 , Reza 18/7/2006

File size: 83.5 KB
Line 
1\documentclass[twoside,10pt]{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
20\usepackage[ps2pdf,bookmarks,bookmarksnumbered,%
21 urlcolor=blue,citecolor=blue,linkcolor=blue,%
22 pagecolor=blue,%hyperindex,%
23 colorlinks=true,hyperfigures=true,hyperindex=true
24 ]{hyperref}
25
26\makeindex % Constitution d'index
27
28\newcommand{\rond}{$\bullet \ $}
29\newcommand{\etoile}{$\star \ $}
30\newcommand{\cercle}{$\circ \ $}
31\newcommand{\carre}{$\Box \ $}
32
33
34\begin{document}
35
36\begin{titlepage}
37% The title page - top of the page with the title of the paper
38\titrehp{Sophya \\ An overview }
39% Authors list
40\auteurs{
41R. Ansari & ansari@lal.in2p3.fr \\
42E. Aubourg & aubourg@hep.saclay.cea.fr \\
43G. Le Meur & lemeur@lal.in2p3.fr \\
44C. Magneville & cmv@hep.saclay.cea.fr \\
45S. Henrot-Versille & versille@in2p3.fr
46}
47% \auteursall
48% The title page - bottom of the page with the paper number
49\vspace{1cm}
50\begin{center}
51{\bf \Large Sophya Version: 2.0 (V\_Jul2006) }
52% Document revision 1.0 }
53\end{center}
54\titrebp{1}
55\end{titlepage}
56
57\tableofcontents
58
59\newpage
60
61\section{Introduction}
62
63{\bf SOPHYA} ({\bf SO}ftware for {\bf PHY}sics {\bf A}nalysis)
64is a collection of C++ classes designed for numerical and
65physics analysis software development. Our goal is to provide
66easy to use, yet powerful classes which can be used by scientists.
67Although some of the SOPHYA modules (SkyMap, Samba, SkyT)
68have been designed with the specific goal CMB data analysis, most
69modules presented here have a much broader scope and can be
70used in scientific data analysis and modeling/simulation.
71Whenever possible, we use existing numerical packages and libraries,
72encapsulating them in classes in order to facilitate their usage.
73\par
74Our main requirements in designing SOPHYA classes can be summarized as
75follow:
76\begin{itemize}
77\item[\rond] Provide a comprehensive set of data containers, such as arrays and tables
78(tuple) covering the most common needs in scientific simulation and data analysis
79softwares.
80\item[\rond] Take advantage of the C++ language and define methods and operators
81for most basic operation, such as arithmetic operations, in a rather intuitive way, while
82maintaining performances comparable to low level coding in other languages
83(C, Fortran, F90 \ldots)
84\item[\rond] Simplify memory management for programmers using the class library.
85This has been a strong requirement for most SOPHYA classes. Automatic reference
86sharing and memory management is implemented in SOPHYA classes intended
87for large size objects. We recommend to allocate SOPHYA objects on the stack,
88including when objects are returned by methods or functions.
89See section \ref{memgt} for more information.
90\item[\rond] Archiving, importing (reading) and exporting (writing) data in a
91efficient and consistent way is a major concern in many scientific software
92and projects. SOPHYA provide a native data I/O or persistence system,
93(PPF, \ref{ppfdesc}) as well as import/export services for ASCII and FITS formats.
94\end{itemize}
95
96% \vspace*{2mm}
97This documents
98presents only a brief overview of the class library,
99mainly from the user's point of view. A more complete description
100can be found in the reference manual, available from the SOPHYA
101web site: % {\bf http://www.sophya.org}.
102\href{http://www.sophya.org}{http://www.sophya.org}.
103%%%
104%%%
105\subsection{SOPHYA modules}
106The source directory tree
107\footnote{ CVS: cvsserver.lal.in2p3.fr:/exp/eros/CVSSophya}
108is organised into a number of modules.
109
110\begin{itemize}
111\item[] {\bf BuildMgr/} Scripts for code management,
112makefile generation and software installation
113\item[] {\bf BaseTools/} General architecture support classes such
114as {\tt PPersist, NDataBlock<T>}, and few utility classes
115such as the dynamic variable list manager ({\tt DVList}) as well
116as the basic set of exception classes used in SOPHYA.
117\item[] {\bf TArray/} template numerical arrays, vectors and matrices \\
118({\tt TArray<T> TMatrix<T> TVector<T> } \ldots)
119\item[] {\bf NTools/} Some standard numerical analysis tools
120(linear, and non linear parameter fitting, FFT, \ldots)
121\item[] {\bf HiStats/} Histogram-ming and data set handling classes (tuples) \\
122({\tt Histo Histo2D NTuple DataTable} \ldots)
123\item[] {\bf SkyMap/} Local and full sky maps, and some 3D geometry
124handling utility classes. \\
125({\tt PixelMap<T>, LocalMap<T>, SphericalMap<T>, \ldots})
126\item[] {\bf SUtils/} This module contains few utility classes, such as the
127{\tt DataCard} class, as well as string manipulation functions in C and C++.
128\item[] {\bf SysTools/} This module contains classes implementing
129an interface to various OS specific services, such
130threads and dynamic link/shared library handling.
131
132\end{itemize}
133
134The modules listed below are more tightly related to the
135CMB (Cosmic Microwave Background) data analysis problem:
136\begin{itemize}
137\item[] {\bf SkyT/}
138classes for spectral emission and detector frequency response modelling \\
139({\tt SpectralResponse, RadSpectra, BlackBody} \ldots)
140\item[] {\bf Samba/} Spherical harmonic analysis, noise generators \ldots
141\end{itemize}
142
143The following modules contain the interface classes with
144external libraries:
145\begin{itemize}
146\item[] {\bf FitsIOServer/} Classes for handling file input-output
147in FITS format using the cfitsio library.
148FITS is maintained by NASA and SAO and is available from: \\
149\href{http://heasarc.gsfc.nasa.gov/docs/software/fitsio/fitsio.html}
150{http://heasarc.gsfc.nasa.gov/docs/software/fitsio/fitsio.html}
151\item[] {\bf LinAlg/} Interface with Lapack linear algebra package.
152Lapack is a linear algebra package and can be downloaded from: \\
153\href{http://www.netlib.org/lapack/}{http://www.netlib.org/lapack/}
154\item[] {\bf IFFTW/} Interface with FFTW package (libfftw.a).
155FFTW is a package for performing Fourier transforms, written in C.
156Documentation and source code can be found at: \\
157\href{http://www.fftw.org/}{http://www.fftw.org/}
158\item[] {\bf XAstroPack/} Interface to some common astronomical
159computation libraries. Presently, this module uses an external library
160extracted from the {\bf Xephem } source code. The corresponding source
161code is also available from SOPHYA cvs repository, module {\bf XephemAstroLib}.
162Information on Xephem can be found at : \\
163\href{http://www.clearskyinstitute.com/xephem/}{http://www.clearskyinstitute.com/xephem/}
164\item[] {\bf MinuitAdapt/} Wrapper classes to CERN minimization routines (Minuit). \\
165\href{http://wwwinfo.cern.ch/asdoc/minuit/minmain.html}{http://wwwinfo.cern.ch/asdoc/minuit/minmain.html}
166
167\end{itemize}
168
169The following modules contain each a set of related programs using the
170SOPHYA library.
171\begin{itemize}
172\item[] {\bf Tests/} Simple test programs
173\item[] {\bf PrgUtil/} Various utility programs (runcxx, scanppf, scanfits, \ldots)
174\item[] {\bf PrgMap/} Programs performing operations on skymaps: projections,
175power spectrum in harmonic space, \ldots
176\item[] {\bf PMixer/} skymixer and related programs
177\end{itemize}
178
179As a companion to SOPHYA, the {\bf (s)piapp} interactive data analysis
180program is built on top of SOPHYA and the {\bf PI} GUI class library
181and application framework. The {\bf PI} ({\bf P}eida {\bf Interactive})
182development started in 1995, in the EROS \footnote{EROS: {\bf E}xp\'erience
183de {\bf R}echerche d'{\bf O}bjets {\bf S}ombres - http://eros.in2p3.fr}
184microlensing search collaboration, with PEIDA++ \footnote {PEIDA++:
185The EROS data analysis class library -
186http://www.lal.in2p3.fr/recherche/eros/PeidaDoc/}.
187The {\bf PI} documentation and the {\bf piapp} user's guide are available
188from \href{http://www.sophya.org}{http://www.sophya.org}.
189%\href{http://www.sophya.org}{http://www.sophya.org}.
190The {\bf PI} is organized as the following modules:
191\begin{itemize}
192\item[] {\bf PI/} Portable GUI class library and application development
193framework kernel.
194\item[] {\bf PIGcont/} Contour-plot drawing classes.
195\item[] {\bf PIext/} Specific drawers and adapters for SOPHYA objects,
196and the {\bf piapp} interactive data analysis framework.
197\item[] {\bf ProgPI/} interactive analysis tool main program and pre-loaded
198modules.
199\end{itemize}
200
201Modules containing examples and demo programs and scripts:
202\begin{itemize}
203\item[] {\bf Examples/} Sample SOPHYA codes and example programs and
204makefiles.
205\item[] {\bf DemoPIApp/} Sample scripts and programs for (s)piapp
206interactive analysis tools.
207\end{itemize}
208\newpage
209
210\section{Using Sophya}
211The organisation of SOPHYA directories and some of the associated
212utility programs are described in this section.
213Basic usage of Sophya classes are described in in the following sections.
214Complete Sophya documentation can be found at our web site
215{\bf http://www.sophya.org}.
216
217\subsection{Directories, environment variables, configuration files}
218\label{directories}
219The environment variable {\bf SOPHYABASE} is used
220to define the path where the Sophya libraries and binaries are installed.
221\begin{itemize}
222\item \$SOPHYABASE/include : Include (.h) files
223\item \$SOPHYABASE/lib : Path for the archive libraries (.a)
224\item \$SOPHYABASE/slb: Shared library path (.so or .dylib on Darwin/MacOS)
225\item \$SOPHYABASE/exe : Path for binary program files
226\end{itemize}
227
228The directory { \tt \$SOPHYABASE/include/SophyaConfInfo/ } contains files
229describing the installed configuration of SOPHYA software.
230
231The file { \tt \$SOPHYABASE/include/machdefs.h } contains definitions
232(flags, typedef) used in SOPHYA, while some more specific flags,
233are found in { \tt \$SOPHYABASE/include/sspvflags.h }
234
235The file { \tt \$SOPHYABASE/include/sophyamake.inc } contains the
236compilation commands and flags used for building the software.
237Users can use most of compilation and link commands defined in this file:
238 {\tt \$CCOMPILE , \$CXXCOMPILE . \$CXXLINK \ldots}.
239 (See module Example).
240
241The configure script (BuildMgr/configure) creates the directory tree and the
242above files. It also copy (or create symbolic link) for all SOPHYA include
243files as well as symbolic links for external libraries
244include files path in {\tt \$SOPHYABASE/include} (FitsIO, FFTW, XAstro \ldots).
245
246Object files for each module are grouped in a static archive library
247by the build procedure (libXXX.a for module
248XXX, with XXX = BaseTools, TArray, HiStats, FitsIOServer \ldots).
249
250When shared libraries are build, all stand alone SOPHYA modules
251are grouped in {\tt libsophya.so}, {\tt libextsophya.so} contains
252the interface modules with external libraries {\bf (FitsIOServer, LinAlg \ldots)},
253while {\bf PI, PIext, PIGcont} modules are grouped in {\tt libPI.so}.
254Alternatively, it is possible to group all modules in a single shared
255library {\tt libAsophyaextPI.so} (See \ref{build})
256
257In order to use the shared libraries, the {\bf LD\_LIBRARY\_PATH} variable
258should contain the Sophya shared library path
259({\tt \$SOPHYABASE/slb}).
260On Silicon Graphics machines with IRIX64 operating system,
261the default SOPHYA configuration correspond to the 64 bit architecture.
262The environment variable { \bf LD\_LIBRARY64\_PATH } replace in
263this case the usual {\bf LD\_LIBRARY\_PATH} variable.
264On IBM machines with AIX, the {\bf LIBPATH} environment variables
265contains the shared libraries search path.
266
267When using the dynamic load services in SOPHYA ({\tt PDynLinkMgr}
268class), in runcxx or (s)piapp applications for example, the shared
269library search path must contain the current working directory (
270dot . in unix).
271
272\subsection{the runcxx program}
273\index{runcxx} \label{runcxx}
274{\bf runcxx} is a simple program which can be used to compile, link
275and run simple C++ programs. It handles the creation of a
276complete program file, containing the basic set C++ include files,
277the necessary include files for SOPHYA SysTools, TArray, HiStats
278and NTools modules, and the main program with exception handling.
279Other Sophya modules can be included using the {\tt -import} flag.
280Use of additional include files can be specified using the
281{\tt -inc} flag.
282\begin{verbatim}
283csh> runcxx -h
284 PIOPersist::Initialize() Starting Sophya Persistence management service
285SOPHYA Version 1.9 Revision 0 (V_Mai2005) -- May 31 2005 15:11:32 cxx
286 runcxx : compiling and running of a piece of C++ code
287 Usage: runcxx [-compopt CompileOptions] [-linkopt LinkOptions]
288 [-tmpdir TmpDirectory] [-f C++CodeFileName]
289 [-inc includefile] [-inc includefile ...]
290 [-import modulename] [-import modulename ...]
291 [-uarg UserArg1 UserArg2 ...]
292 if no file name is specified, read from standard input
293 modulenames: SkyMap, Samba, SkyT, FitsIOServer,
294 LinAlg, IFFTW, XAstroPack
295\end{verbatim}
296Most examples in this manual can be tested using runcxx. The
297example below shows how to compile, link and run a sample
298code.
299\begin{verbatim}
300// File example.icc
301Matrix a(3,3);
302a = IdentityMatrix(1.);
303cout << a ;
304// Executing this sample code
305csh> runcxx -f example.icc
306\end{verbatim}
307
308\subsection{the scanppf program}
309\index{scanppf} \label{scanppf}
310{\bf scanppf} is a simple SOPHYA application which can be used to check
311PPF files and list their contents. It can also provide the list of all registered
312PPF handlers.
313\begin{verbatim}
314csh> scanppf -h
315 PIOPersist::Initialize() Starting Sophya Persistence management service
316SOPHYA Version 2.0 Revision 0 (V_Jul2006) -- Jul 17 2006 14:13:27 cxx
317 Usage: scanppf [flags] filename
318 flags = -s -n -a0 -a1 -a2 -a3 -lh -lho -lmod
319 -s[=default} : Sequential reading of objects
320 -n : Object reading at NameTags
321 -a0...a3 : Tag List with PInPersist.AnalyseTags(0...3)
322 -lh : List PPersist handler classes
323 -lho : List PPersist handler and dataobj classes
324 -lmod : List initialized/registered modules
325\end{verbatim}
326
327\subsection{the scanfits program}
328\index{scanfits} \label{scanfits}
329{\bf scanfits} is a SOPHYA program using the FitsIOServer
330\footnote{FitsIOServer module uses the cfitsio library. scanfits has to be linked with
331with FitsIOServer module and cfitsio libraries, or libextsophya.so}
332module which can be used
333to analyse the content of FITS files. It can list the FITS headers, the appropriate
334SOPHYA-FITS handler (implementing {\tt FitsHandlerInterface}) class, and the list of
335all registered FITS handlers.
336\begin{verbatim}
337csh> scanfits -h
338 PIOPersist::Initialize() Starting Sophya Persistence management service
339SOPHYA Version 2.0 Revision 0 (V_Jul2006) -- Jul 17 2006 14:13:27 cxx
340 Usage: scanfits [flags] filename
341 flags = -V1 -lh -rd -header
342 -V1 : Scan using old (V1) code version
343 -lh : Print the list of registered handlers (FitsHandlerInterface)
344 -rd : try to read each HDU data using appropriate handler
345 -header : List header information
346\end{verbatim}
347
348\newpage
349
350\section{Copy constructor and assignment operator}
351\label{memgt}
352In C++, objects can be copied by assignment or by initialisation.
353Copying by initialisation corresponds to creating an object and
354initialising its value through the copy constructor.
355The copy constructor has its first argument as a reference, or
356const reference to the object's class type. It can have
357more arguments, if default values are provided.
358Copying by assignment applies to an existing object and
359is performed through the assignment operator (=).
360The copy constructor implements this for identical type objects:
361\begin{verbatim}
362class MyObject {
363public:
364 MyObject(); // Default constructor
365 MyObject(MyObject const & a); // Copy constructor
366 MyObject & operator = (MyObject const & a) // Assignment operator
367}
368\end{verbatim}
369The copy constructors play an important role, as they are
370called when class objects are passed by value,
371returned by value, or thrown as an exception.
372\begin{verbatim}
373// A function declaration with an argument of type MyObject,
374// passed by value, and returning a MyObject
375MyObject f(MyObject x)
376{
377 MyObject r;
378 ...
379 return(r); // Copy constructor is called here
380}
381// Calling the function :
382MyObject a;
383f(a); // Copy constructor called for a
384\end{verbatim}
385It should be noted that the C++ syntax is ambiguous for the
386assignment operator. {\tt MyObject x; x=y; } and
387{\tt MyObject x=y;} have different meaning.
388\begin{verbatim}
389MyObject a; // default constructor call
390MyObject b(a); // copy constructor call
391MyObject bb = a; // identical to bb(a) : copy constructor call
392MyObject c; // default constructor call
393c = a; // assignment operator call
394\end{verbatim}
395
396As a general rule in SOPHYA, objects which implements
397reference sharing on their data members have a copy constructor
398which shares the data, while the assignment operator copies or
399duplicate the data.
400
401\newpage
402\section{Module BaseTools}
403
404{\bf BaseTools} contains utility classes such as
405{\tt DVlist}, an hierarchy of exception classes for Sophya, a template
406class {\tcls{NDataBlock}} for handling reference counting on numerical
407arrays, as well as classes providing the services for implementing simple
408serialisation.
409\vspace*{5mm}
410
411\subsection{Initialisation}
412\index{SophyaInitiator}
413A number of actions have to be taken before
414some of the services provided by SOPHYA become operational. This is the case
415of SOPHYA persistence, as well as FITS I/O facilities.
416Initialisation of many SOPHYA modules is performed through an initialiser class,
417which inherits from {\bf SophyaInitiator}. Static instance of each
418initialiser class exist in the library and the various SOPHYA services
419should be operational when the user code ({\tt main()}) starts.
420However, in some cases the run time loader may not perform correctly the static
421object initialisation. You may thus need to instanciate the initialiser class
422for the modules you use in the beginning of your main program: \\
423{\tt TArrayInitiator , HiStatsInitiator , SkyMapInitiator \ldots }
424%%%
425\subsection{SOPHYA persistence}
426\label{ppfdesc}
427\index{PPersist} \index{PInPersist} \index{POutPersist}
428\begin{figure}[hbt]
429\dclsa{PPersist}
430\dclsccc{PPFBinarIOStream}{PPFBinaryInputStream}{PInPersist}
431\dclscc{PPFBinaryOutputStream}{POutPersist}
432\caption{partial class diagram for classes handling persistence in Sophya}
433\end{figure}
434A simple persistence mechanism is defined in SOPHYA. Its main
435features are:
436\begin{itemize}
437\item[] Portable file format, containing the description of the data structures
438and object hierarchy. \\
439{\bf PPF} {\bf P}ortable {\bf P}ersistence file {\bf F}ormat.
440\index{PPF}
441\item[] Handling of read/write for multiply referenced objects.
442\item[] All write operations are carried using sequential access only. This
443holds also for read operations, unless positional tags are used.
444SOPHYA persistence services can thus be used to transfer objects
445through network links.
446\item[] The serialisation (reading/writing) for objects for a given class
447is implemented through a handler object. The handler class inherits
448from {\tt PPersist} class.
449\item[] A run time registration mechanism is used in conjunction with
450RTTI (Run Time Type Identification) for identifying handler classes
451when reading {\bf PInPersist} streams, or for associating handlers
452with data objects {\bf AnyDataObject} for write operations.
453\end{itemize}
454A complete description of SOPHYA persistence mechanism and guidelines
455for writing delegate classes for handling object persistence is beyond
456the scope of this document. The most useful methods for using Sophya
457persistence are listed below:
458\begin{itemize}
459\item[] {\tt POutPersist::PutObject(AnyDataObj \& o)} \\
460Writes the data object {\bf o} to the output stream.
461\item[] {\tt POutPersist::PutObject(AnyDataObj \& o, string tagname)} \\
462Writes the data object {\bf o} to the output stream, associated with an
463identification tag {\bf tagname}.
464\item[] {\tt PInPersist::GetObject(AnyDataObj \& o)} \\
465Reads the next object in stream into {\bf o}. An exception is
466generated for incompatible object types.
467\item[] {\tt PInPersist::GetObject(AnyDataObj \& o, string tagname)} \\
468Reads the object associated with the tag {\bf tagname} into {\bf o}.
469An exception is generated for incompatible object types.
470\end{itemize}
471The operators {\tt operator << (POutPersist ...) } and
472{\tt operator >> (PInPersist ...) } are often overloaded
473to perform {\tt PutObject()} and {\tt GetObject()} operations.
474the {\bf PPFNameTag} (ppfnametag.h) class can be used in conjunction with
475{\tt << >> } operators to write objects with a name tag or to retrieve
476an object identified with a name tag. The example below shows the
477usage of these operators:
478\begin{verbatim}
479// Creating and filling a histogram
480Histo hw(0.,10.,100);
481...
482// Writing histogram to a PPF stream
483POutPersist os("hw.ppf");
484os << PPFNameTag("myhisto") << hw;
485
486// Reading a histogram from a PPF stream
487PInPersist is("hr.ppf");
488is >> PPFNameTag("myhisto") >> hr;
489\end{verbatim}
490
491The {\bf scanppf} program can be used to list the content of a PPF file.
492
493\subsection{\tcls{NDataBlock}}
494\index{\tcls{NDataBlock}}
495\begin{figure}[hbt]
496\dclsbb{AnyDataObj}{\tcls{NDataBlock}}
497\dclsbb{PPersist}{\tcls{FIO\_NDataBlock}}
498\end{figure}
499The {\bf \tcls{NDataBlock}} is designed to handle reference counting
500and sharing of memory blocs (contiguous arrays) for numerical data
501types. Initialisation, resizing, basic arithmetic operations, as
502well as persistence handling services are provided.
503The persistence handler class ({\tt \tcls{FIO\_NDataBlock}}) insures
504that a single copy of data is written for multiply referenced objects,
505and the data is shared among objects when reading.
506\par
507The example below shows writing of NDataBlock objects through the
508use of overloaded operator $ << $ :
509\begin{verbatim}
510#include "fiondblock.h"
511// ...
512POutPersist pos("aa.ppf");
513NDataBlock<r_4> rdb(40);
514rdb = 567.89;
515pos << rdb;
516// We can also use the PutObject method
517NDataBlock<int_4> idb(20);
518idb = 123;
519pos.PutObject(idb);
520\end{verbatim}
521The following sample programs show the reading of the created PPF file :
522\begin{verbatim}
523PInPersist pis("aa.ppf");
524NDataBlock<r_4> rdb;
525pis >> rdb;
526cout << rdb;
527NDataBlock<int_4> idb;
528cout << idb;
529\end{verbatim}
530
531\subsection{DVList, MuTyV and TimeStamp classes}
532\index{DVList} \index{MuTyV} \index{TimeStamp}
533\begin{figure}[hbt]
534\dclsa{MuTyV}
535\dclsbb{AnyDataObj}{DVList}
536\dclsbb{PPersist}{\tclsc{ObjFileIO}{DVList}}
537\end{figure}
538The {\bf DVList} class objects can be used to create and manage list
539of values, associated with names. A list of pairs of (MuTyV, name(string))
540is maintained by DVList objects. {\bf MuTyV} is a simple class
541capable of holding string, integer, float or complex values,
542providing easy conversion methods between these objects.
543{\bf MuTyV} objects can also hold {\bf TimeStamp } objects.
544\begin{verbatim}
545// Using MuTyV objects
546MuTyV s("hello"); // string type value
547MuTyV x;
548x = "3.14159626"; // string type value, ASCII representation for Pi
549double d = x; // x converted to double = 3.141596
550x = 314; // x contains the integer value = 314
551// Using DVList
552DVList dvl;
553dvl("Pi") = 3.14159626; // float value, named Pi
554dvl("Log2") = 0.30102999; // float value, named Log2
555dvl("FileName") = "myfile.fits"; // string value, named myfile.fits
556// Printing DVList object
557cout << dvl;
558\end{verbatim}
559
560\begin{figure}[hbt]
561\dclsbb{AnyDataObj}{TimeStamp}
562\end{figure}
563%
564The {\bf TimeStamp} class represent date and time and provides
565many standard operations, such as Initialisation from strings,
566conversion to strings and time interval computations. \\
567Usage example:
568\begin{verbatim}
569// Create a object with the current date and time and prints it to cout
570TimeStamp ts;
571cout << ts << endl;
572// Create an object with a specified date and time
573TimeStamp ts2("01/01/1905","00:00:00");
574// Get the number of days since 0 Jan 1901
575cout << ts2.ToDays() << endl;
576
577// Combined use of TimeStamp and MuTyV
578string s;
579TimeStamp ts; // Current date/time
580MuTyV mvt = ts;
581s = mvt; // s contains the current date in string format
582cout << s << endl;
583\end{verbatim}
584
585\subsection{\tcls{SegDataBlock} , \tcls{SwSegDataBlock}}
586\index{\tcls{SegDataBlock}} \index{\tcls{SwSegDataBlock}}
587%%
588\begin{figure}[hbt]
589\dclsccc{AnyDataObj}{\tcls{SegDBInterface}}{ \tcls{SegDataBlock} }
590\dclscc{\tcls{SegDBInterface}}{ \tcls{SwSegDataBlock} }
591\end{figure}
592\begin{itemize}
593\item[] \tcls{SegDataBlock} handles arrays of object of
594type {\bf T} with reference sharing in memory. The array can be extended
595(increase in array size) with fixed segment size. It implements the interface
596defined by \tcls{SegDBInterface}.
597\item[] \tcls{SwSegDataBlock} Implements the same \tcls{SegDBInterface}
598using a data swapper object. Data swappers implement the interface defined in
599(\tcls{DataSwapperInterface} class. \tcls{SwSegDataBlock} can
600thus be used to handle arrays with very large number of objects.
601These classes handles reference sharing.
602\end{itemize}
603
604\newpage
605\section{Module TArray}
606\index{\tcls{TArray}}
607{\bf TArray} module contains template classes for handling standard
608operations on numerical arrays. Using the class {\tt \tcls{TArray} },
609it is possible to create and manipulate up to 5-dimension numerical
610arrays {\tt (int, float, double, complex, \ldots)}. The include
611file {\tt array.h} declares all the classes and definitions
612in module TArray. {\bf Array} is a typedef for arrays
613with double precision floating value elements. \\
614{\tt typedef TArray$<$r\_8$>$ Array ; }
615
616\begin{figure}[hbt]
617\dclsccc{AnyDataObj}{BaseArray}{\tcls{TArray}}
618\dclsbb{PPersist}{\tcls{FIO\_TArray}}
619\end{figure}
620
621The development of this module started around 1999-2000,
622after evaluation of a number of publicly available
623C++ array hadling packages, including TNT, Lapack++, Blitz++,
624as well as commercial packages from RogueWave (math.h++ \ldots).
625Most of these packages provide interesting functionalities, however,
626not any one package seemed to fulfill most of our requirements.
627\begin{itemize}
628\item Capability to handle {\bf large - multidimensional - dense}
629arrays, for numerical data types. Although we have used templates, for
630data type specialisation, the actual code, apart inline functions is
631not in header files. Instead, we use explicit instanciation, and the
632compiled code for the various numerical types of arrays is the
633library .
634\item The shape and size of the arrays can be defined and changed
635at run time. The classes ensure the memory management of the
636created objets, with reference sharing for the array data.
637The default behaviour of the copy constructor is to share the data,
638avoiding expensive memory copies.
639\item The package provides transparent management of sub-arrays
640and slices, in an intuitive way, somehow similar to what is
641available in Mathlab or Scilab.
642\item The memory organisation for arrays, specially matrices
643(row-major or column major) can be
644controled. This provide compatibility when using existing C or
645Fortran coded numerical libraries.
646\item The classes provide efficient methods to perform basic arithmetic
647and mathematical operations on arrays. In addition, operator overload
648provides intuitive programming for element acces and most basic
649arithmetic operations.
650\item Conversion can be performed between arrays with different
651data types. Copy and arithmetic operations can be done transparently
652between arrays with different memory organisation patterns.
653\item This module does not provide more complex operations
654such as FFT or linear algebra. Additional libraries are used, with interface
655classes for these operations.
656\item ASCII formatted I/O, for printing and read/write operations to/from text files.
657\item Efficient binary I/O for object persistence (PPF format), or import/export
658to other data formats, such as FITS are provided by helper or handler classes.
659\end{itemize}
660
661\subsection{Using arrays}
662\index{Sequence} \index{RandomSequence} \index{RegularSequence}
663\index{EnumeratedSequence}
664The example below shows basic usage of arrays, creation, initialisation
665and arithmetic operations. Different kind of {\bf Sequence} objects
666can be used for initialising arrays.
667
668\begin{figure}[hbt]
669\dclsbb{Sequence}{RandomSequence}
670\dclsb{RegularSequence}
671\dclsb{EnumeratedSequence}
672\end{figure}
673
674The example below shows basic usage of arrays:
675\index{\tcls{TArray}}
676\begin{verbatim}
677// Creating and initialising a 1-D array of integers
678TArray<int> ia(5);
679EnumeratedSequence es;
680es = 24, 35, 46, 57, 68;
681ia = es;
682cout << "Array<int> ia = " << ia;
683// 2-D array of floats
684TArray<r_4> b(6,4), c(6,4);
685// Initializing b with a constant
686b = 2.71828;
687// Filling c with random numbers
688c = RandomSequence();
689// Arithmetic operations
690TArray<r_4> d = b+0.3f*c;
691cout << "Array<float> d = " << d;
692\end{verbatim}
693
694The copy constructor shares the array data, while the assignment operator
695copies the array elements, as illustrated in the following example:
696\begin{verbatim}
697TArray<int> a1(4,3);
698a1 = RegularSequence(0,2);
699// Array a2 and a1 shares their data
700TArray<int> a2(a1);
701// a3 and a1 have the same size and identical elements
702TArray<int> a3;
703a3 = a1;
704// Changing one of the a2 elements
705a2(1,1,0) = 555;
706// a1(1,1) is also changed to 555, but not a3(1,1)
707cout << "Array<int> a1 = " << a1;
708cout << "Array<int> a3 = " << a3;
709\end{verbatim}
710
711\subsection{Arithmetic operations}
712The four usual arithmetic operators ({\bf + \, - \, * \, / }) are defined
713to perform constant addition, subtraction, multiplication and division.
714The three operators ({\bf + \, - \, / }) between two arrays of the same type
715are defined to perform element by element addition, subtraction
716and division. In order to avoid confusion with matrix multiplication,
717element by element multiplication is defined by overloading the
718operator {\bf \, \&\& \, }, as shown in the example below:
719\begin{verbatim}
720TArray<int_4> a(4,3), b(4,3), c , d, e;
721a = RegularSequence(1.,1.);
722b = RegularSequence(10.,10.);
723cout << a << b ;
724c = a && b;
725d = c / a;
726e = (c / b) - a;
727cout << c << d << e;
728\end{verbatim}
729
730\subsection{Matrices and vectors}
731\index{\tcls{TMatrix}} \index{\tcls{TVector}}
732\begin{figure}[hbt]
733\dclsccc{\tcls{TArray}}{\tcls{TMatrix}}{\tcls{TVector}}
734\end{figure}
735Vectors and matrices are 2 dimensional arrays. The array size
736along one dimension is equal 1 for vectors. Column vectors
737have {\tt NCols() = 1} and row vectors have {\tt NRows() = 1}.
738Mathematical expressions involving matrices and vectors can easily
739be translated into C++ code using {\tt TMatrix} and
740{\tt TVector} objects. {\bf Matrix} and {\bf Vector} are
741typedefs for double precision float matrices and vectors.
742The operator {\bf *} beteween matrices is redefined to
743perform matrix multiplication. One can then write: \\
744\begin{verbatim}
745 // We create a row vector
746 Vector v(1000, BaseArray::RowVector);
747 // Initialize values with a random sequence
748 v = RandomSequence();
749 // Compute the vector length (norm)
750 double norm = (v*v.Transpose()).toScalar();
751 cout << "Norm(v) = " << norm << endl;
752\end{verbatim}
753
754This module contains basic array and matrix operations
755such as the Gauss matrix inversion algorithm
756which can be used to solve linear systems, as illustrated by the
757example below:
758\begin{verbatim}
759#include "sopemtx.h"
760// ...
761// Creation of a random 5x5 matrix
762Matrix A(5,5);
763A = RandomSequence(RandomSequence::Flat);
764Vector X0(5);
765X0 = RandomSequence(RandomSequence::Gaussian);
766// Computing B = A*X0
767Vector B = A*X0;
768// Solving the system A*X = B
769Vector X;
770LinSolve(A, B, X);
771// Checking the result
772Vector diff = X-X0;
773cout << "X-X0= " << diff ;
774double min,max;
775diff.MinMax(min, max);
776cout << " Min(X-X0) = " << min << " Max(X-X0) = " << max << endl;
777\end{verbatim}
778
779{\bf Warning: } The operations defined in {\tt sopemtx.h}, such as
780matrix inversion and linear system solver use a basic Gauss pivot
781algorithm which are not adapted for large matrices ($>\sim 100x100$).
782The services provided in other modules, such as {\bf LinAlg} should
783be preferred in such cases.
784
785\subsection{Working with sub-arrays and Ranges}
786\index{Range}
787A powerful mechanism is included in array classes for working with
788sub-arrays. The class {\bf Range} can be used to specify range of array
789indexes in any of the array dimensions. Any regularly spaced index
790range can be specified, using the {\tt start} and {\tt end} index
791and an optional step (or stride). It is also possible to specify
792the {\tt start} index and the number of elements.
793\begin{itemize}
794\item {\bf Range::all()} {\tt = Range(Range::firstIndex(), Range::lastIndex())} \\
795return a Range objects representing all valid indexes along the
796corresponding axe.
797\item {\bf Range::first()} {\tt = Range(Range::firstIndex())} \\
798return a Range object representing the first valid index
799\item {\bf Range::last()} {\tt = Range(Range::lastIndex())}
800return a Range object representing the last valid index
801\item {\bf Range(idx)} represents a single index ({\bf = idx})
802\item {\bf Range(first, last)} represents the range of indices
803{\bf first} $\leq$ index $\leq$ {\bf last}.
804The static method {\tt Range::lastIndex()} can be used
805to specify the last valid index.
806\item {\bf Range(first, last, step)} represents the range of index
807which is equivalent to \\ {\tt for(index=first; index <= last; index += step) }
808\item { \bf Range (first, last, size, step) } the general form can be used
809to specify an index range, using the number of elements.
810It is possible to specify a range of index, ending with the last valid index.
811For example \\
812\hspace*{5mm}
813{\tt Range(Range::lastIndex(), Range::lastIndex(), 3, 2) } \\
814defines the index range: \hspace*{5mm} last-4, last-2, last.
815
816\begin{center}
817\begin{tabular}{ll}
818\hline \\
819\multicolumn{2}{c}{ {\bf Range} {\tt (start, end, size, step) } } \\[2mm]
820\hline \\
821{\bf Range} {\tt r(7); } & index range: \hspace{2mm} 7 \\
822{\bf Range} {\tt r(3,6); } & index range: \hspace{2mm} 3,4,5,6 \\
823{\bf Range} {\tt r(3,7,2); } & index range: \hspace{2mm} 3,5,7 \\
824{\bf Range} {\tt r(7,0,3,1); } & index range: \hspace{2mm} 7,8,9 \\
825{\bf Range} {\tt r(10,0,5,2); } & index range: \hspace{2mm} 10,12,14,16,18 \\[2mm]
826\hline
827\end{tabular}
828\end{center}
829\end{itemize}
830
831The method {\tt TArray<T>SubArray(Range ...)} can be used
832to extract subarrays and slices. The operator {\tt operator() (Range rx, Range ry ...)}
833is also overloaded for sub-array extraction.
834For matrices, {\tt TMatrix<T>::Row()} and {\tt TMatrix<T>::Column()}
835extract selected matrix rows and columns.
836
837The example illustrates the sub-array extraction using Range objects:
838\begin{verbatim}
839 // Creating and initialising a 2-D (6 x 4) array of integers
840 TArray<int> iaa(6, 4);
841 iaa = RegularSequence(1,2);
842 cout << "Array<int> iaa = \n" << iaa;
843 // We extract a sub-array - data is shared with iaa
844 TArray<int> iae = iaa(Range(1, Range::lastIndex(), 3) ,
845 Range::all(), Range::first() );
846 cout << "Array<int> iae=subarray(iaa) = \n" << iae;
847 // Changing iae elements changes corresponding iaa elements
848 iae = 0;
849 cout << "Array<int> iae=0 --> iaa = \n" << iaa;
850
851\end{verbatim}
852
853In the following example, a simple low-pass filter, on a one
854dimensional stream (Vector) has been written using sub-arrays:
855
856\begin{verbatim}
857// Input Vector containing a noisy periodic signal
858 Vector in(1024), out(1024);
859 in = RandomSequence(RandomSequence::Gaussian, 0., 1.);
860 for(int kk=0; kk<in.Size(); kk++)
861 in(kk) += 2*sin(kk*0.05);
862// Compute the output vector by a simple low pass filter
863 Vector out(1024);
864 int w = 2;
865 for(int k=w; k<in.Size()-w; k++)
866 out(k) = in(Range(k-w, k+w).Sum()/(2.*w+1.);
867\end{verbatim}
868
869\subsection{Input, Output}
870Arrays can easily be saved to, or restored from files in different formats.
871SOPHYA library can handle array I/O to ASCII formatted files, to PPF streams,
872as well as to files in FITS format.
873FITS format input/output is provided through the classes in
874{\bf FitsIOServer} module. Only arrays with data types
875supported by the FITS standard can be handled during
876I/O operations to and from FITS streams (See the FitsIOServer section
877for additional details).
878
879\subsubsection{PPF streams}
880
881SOPHYA persistence (PPF streams) handles reference sharing, and multiply
882referenced objects are only written once. A hierarchy of arrays and sub-arrays
883written to a PPF stream is thus completely recovered, when the stream is read.
884The following example illustrates this point:
885\begin{verbatim}
886{
887// Saving an array with a sub-array into a POutPersist file
888Matrix A(3,4);
889A = RegularSequence(10,5);
890// Create a sub-array of A
891Matrix AS = A(Range(1,2), Range(2,3));
892// Save the two arrays to a PPF stream
893POutPersist pos("aas.ppf");
894pos << A << AS;
895}
896{
897// Reading arrays from the previously created PPF file aas.ppf
898PInPersist pis("aas.ppf");
899Matrix B,BS;
900pis >> B >> BS;
901// BS is a sub-array of B, modifying BS changes also B
902BS(1,1) = 98765.;
903cout << " B , BS after BS(1,1) = 98765. "
904 << B << BS << endl;
905}
906\end{verbatim}
907The execution of this sample code creates the file {\tt aas.ppf} and
908its output is reproduced here. Notice that the array hierarchy is
909recovered. BS is a sub-array of B, and modifying BS changes also
910the corresponding element in B.
911\begin{verbatim}
912 B , BS after BS(1,1) = 98765.
913
914--- TMatrix<double>(NRows=3, NCols=4) ND=2 SizeX*Y*...= 4x3 ---
91510 15 20 25
91630 35 40 45
91750 55 60 98765
918
919--- TMatrix<double>(NRows=2, NCols=2) ND=2 SizeX*Y*...= 2x2 ---
92040 45
92160 98765
922\end{verbatim}
923
924\centerline{\bf Warning: }
925
926There is a drawback in this behaviour: only a single
927copy of an array is written to a file, even if the array is modified,
928without being resized and written to a PPF stream.
929\begin{verbatim}
930{
931POutPersist pos("mca.ppf");
932TArray<int_4> ia(5,3);
933ia = 8;
934pos << ia;
935ia = 16;
936pos << ia;
937ia = 32;
938pos << ia;
939}
940\end{verbatim}
941
942Only a single copy of the data is effectively written to the output
943PPF file, corresponding to the value 8 for array elements. When we
944read the three array from the file mca.ppf, the same array elements
945are obtained three times (all elements equal to 8):
946\begin{verbatim}
947{
948PInPersist pis("mca.ppf");
949TArray<int_4> ib;
950pis >> ib;
951cout << " First array read from mca.ppf : " << ib;
952pis >> ib;
953cout << " Second array read from mca.ppf : " << ib;
954pis >> ib;
955cout << " Third array read from mca.ppf : " << ib;
956}
957\end{verbatim}
958
959\subsubsection{ASCII streams}
960
961The {\bf WriteASCII} method can be used to dump an array to an ASCII
962formatted file, while the {\bf ReadASCII} method can be used to decode
963ASCII formatted files. Space or tabs are the possible separators.
964Complex numbers should be specified as a pair of comma separated
965real and imaginary parts, enclosed in parenthesis.
966
967\begin{verbatim}
968{
969// Creating array A and writing it to an ASCII file (aaa.txt)
970Array A(4,6);
971A = RegularSequence(0.5, 0.2);
972ofstream ofs("aaa.txt");
973A.WriteASCII(ofs);
974}
975{
976// Decoding the ASCII file aaa.txt
977ifstream ifs("aaa.txt");
978Array B;
979sa_size_t nr, nc;
980B.ReadASCII(ifs,nr,nc);
981cout << " Array B; B.ReadASCII() from file " << endl;
982cout << B ;
983}
984\end{verbatim}
985
986
987\subsection{Complex arrays}
988The {\bf TArray} module provides few functions for manipulating
989arrays of complex numbers (single and double precision).
990These functions are declared in {\tt matharr.h}.
991\begin{itemize}
992\item[\bul] Creating a complex array through the specification of the
993real and imaginary parts.
994\item[\bul] Functions returning arrays corresponding to real and imaginary
995parts of a complex array: {\tt real(za) , imag(za) }
996({\bf Warning:} Note that the present implementation does not provide
997shared memory access to real and imaginary parts.)
998\item[\bul] Functions returning arrays corresponding to the module,
999phase, and module squared of a complex array:
1000 {\tt phase(za) , module(za) , module2(za) }
1001\end{itemize}
1002
1003\begin{verbatim}
1004 TVector<r_4> p_real(10, BaseArray::RowVector);
1005 TVector<r_4> p_imag(10, BaseArray::RowVector);
1006 p_real = RegularSequence(0., 0.5);
1007 p_imag = RegularSequence(0., 0.25);
1008 TVector< complex<r_4> > zvec = ComplexArray(p_real, p_imag);
1009 cout << " :: zvec= " << zvec;
1010 cout << " :: real(zvec) = " << real(zvec) ;
1011 cout << " :::: imag(zvec) = " << imag(zvec) ;
1012 cout << " :::: module2(zvec) = " << module2(zvec) ;
1013 cout << " :::: module(zvec) = " << module(zvec) ;
1014 cout << " :::: phase(zvec) = " << phase(zvec) ;
1015\end{verbatim}
1016
1017The decoding of complex numbers from an ASCII formatted stream
1018is illustrated by the next example. As mentionned already,
1019complex numbers should be specified as a pair of comma separated
1020real and imaginary parts, enclosed in parenthesis.
1021
1022\begin{verbatim}
1023csh> cat zzz.txt
1024(1.,-1) (2., 2.5) -3. 12.
1025-24. (-6.,7.) 14.2 (8.,64.)
1026
1027// Decoding of complex numbers from an ASCII file
1028// Notice that the << operator can be used instead of ReadASCII
1029TArray< complex<r_4> > Z;
1030ifstream ifs("zzz.txt");
1031ifs >> Z;
1032cout << " TArray< complex<r_4> > Z from file zzz.txt " << Z ;
1033\end{verbatim}
1034
1035
1036\subsection{Memory organisation}
1037{\tt \tcls{TArray} } can handle numerical arrays with various memory
1038organisation, as long as the spacing (steps) along each axis is
1039regular. The five axis are labeled X,Y,Z,T,U. The examples below
1040illustrates the memory location for a 2-dimensional, $N_x=4 \times N_y=3$.
1041The first index is along the X axis and the second index along the Y axis.
1042\begin{verbatim}
1043 | (0,0) (0,1) (0,2) (0,3) |
1044 | (1,0) (1,1) (1,2) (1,3) |
1045 | (2,0) (2,1) (2,2) (2,3) |
1046\end{verbatim}
1047In the first case, the array is completely packed
1048($Step_X=1, Step_Y=N_X=4$), with zero offset,
1049while in the second case, $Step_X=2, Step_Y=10, Offset=10$:
1050\begin{verbatim}
1051 | 0 1 2 3 | | 10 12 14 16 |
1052Ex1 | 4 5 6 7 | Ex2 | 20 22 24 26 |
1053 | 8 9 10 11 | | 30 32 34 36 |
1054\end{verbatim}
1055
1056For matrices and vectors, an optional argument ({\tt MemoryMapping})
1057can be used to select the memory mapping, where two basic schemes
1058are available: \\
1059{\tt CMemoryMapping} and {\tt FortranMemoryMapping}. \\
1060In the case where {\tt CMemoryMapping} is used, a given matrix line
1061is packed in memory, while the columns are packed when
1062{\tt FortranMemoryMapping} is used. The first index when addressing
1063the matrix elements (line number index) runs along
1064the Y-axis if {\tt CMemoryMapping} is used, and along the X-axis
1065in the case of {\tt FortranMemoryMapping}.
1066Arithmetic operations between matrices
1067with different memory organisation is allowed as long as
1068the two matrices have the same sizes (Number of rows and columns).
1069The following code example and the corresponding output illustrates
1070these two memory mappings. The {\tt \tcls{TMatrix}::TransposeSelf() }
1071method changes effectively the matrix memory mapping, which is also
1072the case of {\tt \tcls{TMatrix}::Transpose() } method without argument.
1073
1074\begin{verbatim}
1075TArray<r_4> X(4,2);
1076X = RegularSequence(1,1);
1077cout << "Array X= " << X << endl;
1078TMatrix<r_4> X_C(X, true, BaseArray::CMemoryMapping);
1079cout << "Matrix X_C (CMemoryMapping) = " << X_C << endl;
1080TMatrix<r_4> X_F(X, true, BaseArray::FortranMemoryMapping);
1081cout << "Matrix X_F (FortranMemoryMapping) = " << X_F << endl;
1082\end{verbatim}
1083This code would produce the following output (X\_F = Transpose(X\_C)) :
1084\begin{verbatim}
1085Array X=
1086--- TArray<f> ND=2 SizeX*Y*...= 4x2 ---
10871, 2, 3, 4
10885, 6, 7, 8
1089
1090Matrix X_C (CMemoryMapping) =
1091--- TMatrix<f>(NRows=2, NCols=4) ND=2 SizeX*Y*...= 4x2 ---
10921, 2, 3, 4
10935, 6, 7, 8
1094
1095Matrix X_F (FortranMemoryMapping) =
1096--- TMatrix<f>(NRows=4, NCols=2) ND=2 SizeX*Y*...= 4x2 ---
10971, 5
10982, 6
10993, 7
11004, 8
1101\end{verbatim}
1102
1103\newpage
1104
1105\section{Module HiStats}
1106\begin{figure}[hbt]
1107\dclsbb{AnyDataObj}{Histo}
1108\dclscc{Histo}{HProf}
1109\dclscc{Histo}{HistoErr}
1110\dclsbb{AnyDataObj}{Histo2D}
1111\caption{partial class diagram for histograms and ntuples}
1112\end{figure}
1113
1114{\bf HiStats} contains classes for creating, filling, printing and
1115doing various operations on one or two dimensional histograms
1116{\tt Histo} and {\tt Histo2D} as well as profile histograms {\tt HProf}. \\
1117This module also contains {\tt NTuple} and {\tt DataTable} which are
1118more or less the same as the binary or ascii FITS tables.
1119
1120\subsection{Histograms}
1121\subsubsection{1D Histograms}
1122\index{Histo}
1123For 1D histograms, various numerical methods are provided such as
1124computing means and sigmas, finding maxima, fitting, rebinning,
1125integrating \dots \\
1126The example below shows creating and filling a one dimensional histogram
1127of 100 bins from $-5.$ to $+5.$ to create a Gaussian normal distribution
1128with errors~:
1129\begin{verbatim}
1130#include "histos.h"
1131// ...
1132Histo H(-0.5,0.5,100);
1133H.Errors();
1134for(int i=0;i<25000;i++) {
1135 double x = NorRand();
1136 H.Add(x);
1137}
1138H.Print(80);
1139\end{verbatim}
1140
1141\subsubsection{2D Histograms}
1142\index{Histo2D}
1143Much of these operations are also valid for 2D histograms. 1D projection
1144or slices can be set~:
1145\begin{verbatim}
1146#include "histos2.h"
1147// ...
1148Histo2D H2(-1.,1.,100,0.,60.,50);
1149H2.SetProjX(); // create the 1D histo for X projection
1150H2.SetBandX(25.,35.); // create 1D histo projection for 25.<y<35.
1151H2.SetBandX(35.,45.); // create 1D histo projection for 35.<y<45.
1152H2.SetBandX(40.,55.); // create 1D histo projection for 40.<y<55.
1153//... fill H2 with what ever you want
1154H2.Print();
1155Histo *hx = H2.HProjX();
1156 hx->Print(80);
1157Histo *hbx2 = HBandX(1); // Get the second X band (35.<y<45.)
1158 hbx2->Print(80);
1159\end{verbatim}
1160
1161\subsubsection{Profile Histograms}
1162\index{HProf}
1163Profiles histograms {\bf HProf} contains the mean and the
1164sigma of the distribution
1165of the values filled in each bin. The sigma can be changed to
1166the error on the mean. When filled, the profile histogram looks
1167like a 1D histogram and much of the operations that can be done on 1D histo
1168may be applied onto profile histograms.
1169
1170\subsection{Data tables (tuples)}
1171\begin{figure}[hbt]
1172\dclsbb{AnyDataObj}{NTuple}
1173\dclsccc{AnyDataObj}{BaseDataTable}{DataTable}
1174\dclscc{BaseDataTable}{SwPPFDataTable}
1175\end{figure}
1176
1177\subsubsection{NTuple}
1178\index{NTuple}
1179NTuple are memory resident tables of 32 or 64 bits floating values
1180(float/double).They are arranged in columns. Each line is often called an event.
1181These objects are frequently used to analyze data.
1182The piapp graphicals tools can plot a column against an other one
1183with respect to various selection cuts. \\
1184Here is an example of creation and filling~:
1185\begin{verbatim}
1186#include "ntuple.h"
1187#include "srandgen.h"
1188// ...
1189char* nament[4] = {"i","x","y","ey"};
1190r_4 xnt[4];
1191NTuple NT(4,nament);
1192for(i=0;i<5000;i++) {
1193 xnt[0] = i+1;
1194 xnt[1] = 5.*drandpm1(); // a random value between -5 and +5
1195 xnt[2] = 100.*exp(-0.5*xnt[1]*xnt[1]) + 1.;
1196 xnt[3] = sqrt(xnt[2]);
1197 xnt[2] += xnt[3] * NorRand(); // add a random gaussian error
1198 NT.Fill(xnt);
1199}
1200\end{verbatim}
1201
1202XNTuple provide additional functionalities, compared to NTuple.
1203They are deprecated and are only kept for backward compatibility
1204and should not be used anymore. Use DataTable and
1205SwPPFDataTable instead.
1206Object of type XNTuple handle various types
1207of column values (double,float,int,string,...) and can handle
1208very large data sets, through swap space on disk.
1209
1210\subsubsection{DataTables}
1211\label{datatables}
1212\index{DataTable}
1213The class {\bf DataTable} extends significantly the functionalities provided by
1214NTuple. DataTable is a memory resident implementation of the interface
1215{\bf BaseDataTable } which organizes the data as a 2-D table. User can define
1216the name and data type of each column. Data is added to the table as rows.
1217The table is extended as necessary when adding rows.
1218The sample code below shows an example of DataTable usage :
1219\begin{verbatim}
1220 #include "datatable.h"
1221 // ...
1222 {
1223 DataTable dt(64);
1224 dt.AddFloatColumn("X0_f");
1225 dt.AddFloatColumn("X1_f");
1226 dt.AddDoubleColumn("X0X0pX1X1_d");
1227 double x[5];
1228 for(int i=0; i<63; i++) {
1229 x[0] = (i/9)-4.; x[1] = (i/9)-3.; x[2] = x[0]*x[0]+x[1]*x[1];
1230 dt.AddLine(x);
1231 }
1232 // Printing table info
1233 cout << dt ;
1234 // Saving object into a PPF file
1235 POutPersist po("dtable.ppf");
1236 po << dt ;
1237 }
1238\end{verbatim}
1239
1240
1241\index{SwPPFDataTable}
1242The class {\bf SwPPFDataTable} implements the BaseDataTable interface
1243using segmented data blocks with swap on PPF streams. Very large data sets
1244can be created and manipulated through tis class
1245
1246\index{DataTableRow}
1247The class {\bf DataTableRow } is an auxiliary class which simplifies the manipulation
1248of BaseDataTable object rows.
1249The example below show how to create and filling a table, using a PPF stream as
1250swap space. In addition, we have used a {\tt DataTableRow} to prepare data
1251for each table line.
1252\begin{verbatim}
1253 #include "swppfdtable.h"
1254 // ...
1255 {
1256 // ---------- Create an output PPF stream (file)
1257 POutPersist po("swdtable.ppf");
1258 // ------------------
1259 // Create a table with 3 columns, using the above stream as swap space
1260 SwPPFDataTable dtrow(po, 64);
1261 dtrow.AddStringColumn("sline");
1262 dtrow.AddIntegerColumn("line");
1263 dtrow.AddDateTimeColumn("datime");
1264 //
1265 TimeStamp ts, ts2; // Initialize current date and time
1266 string sline;
1267 //---- Create a table row with the required structure
1268 DataTableRow row = dtrow.EmptyRow();
1269 // ----- Fill the table
1270 for(int k = 0; k<2500; k++) {
1271 sline = "L-";
1272 sline += (string)MuTyV(k);
1273 row["sline"] = sline;
1274 row[1] = k;
1275 ts2.Set(ts.ToDays()+(double)k);
1276 row["datime"] = ts2;
1277 dtrow.AddRow(row);
1278 }
1279 //------ Write the table itself to the stream, before closing the file
1280 po << PPFNameTag("SwTable") << dtrow;
1281 }
1282\end{verbatim}
1283%%
1284The previously created table can easily be read in, as shown below:
1285%%
1286\begin{verbatim}
1287 #include "swppfdtable.h"
1288 // ...
1289 {
1290 // ------ Create the input PPF stream (file)
1291 PInPersist pin("swdtable.ppf");
1292 // ------ Read in the SwPPFDataTable object
1293 SwPPFDataTable dtr;
1294 pin >> PPFNameTag("SwTable") >> dtr;
1295 // ---- Create a table row with the required structure
1296 DataTableRow row = dtr.EmptyRow();
1297 // ---- Acces and print two of the table rows :
1298 cout << dtr.GetRow(6, row) << endl;
1299 cout << dtr.GetRow(17, row) << endl;
1300 }
1301\end{verbatim}
1302
1303\subsection{Writing, viewing \dots }
1304
1305All these objects have been design to be written to or read from a persistent file.
1306The following example shows how to write the previously created objects
1307into such a file~:
1308\begin{verbatim}
1309//-- Writing
1310{
1311char *fileout = "myfile.ppf";
1312string tag;
1313POutPersist outppf(fileout);
1314tag = "H"; outppf.PutObject(H,tag);
1315tag = "H2"; outppf.PutObject(H2,tag);
1316tag = "NT"; outppf.PutObject(NT,tag);
1317} // closing ``}'' destroy ``outppf'' and automatically close the file !
1318\end{verbatim}
1319
1320Sophya graphical tools (spiapp) can automatically display and operate
1321all these objects.
1322
1323\newpage
1324\section{Module NTools}
1325
1326This module provides elementary numerical tools for numerical integration,
1327fitting, sorting and ODE solving. FFTs are also provided (Mayer,FFTPack).
1328
1329\subsection{Fitting}
1330\index{Fitting} \index{Minimisation}
1331Fitting is done with two classes {\tt GeneralFit} and {\tt GeneralFitData}
1332and is based on the Levenberg-Marquardt method.
1333\index{GeneralFit} \index{GeneralFitData}
1334GeneralFitData is a class which provide a description of the data
1335to be fitted. GeneralFit is the fitter class. Parametrized functions
1336can be given as classes which inherit {\tt GeneralFunction}
1337or as simple C functions. Classes of pre-defined functions are provided
1338(see files fct1dfit.h and fct2dfit.h). The user interface is very close
1339from that of the CERN {\tt Minuit} fitter.
1340Number of objects (Histo, HProf \dots ) are interfaced with GeneralFit
1341and can be easily fitted. \\
1342Here is a very simple example for fitting the previously created NTuple
1343with a Gaussian~:
1344\begin{verbatim}
1345#include "fct1dfit.h"
1346// ...
1347
1348// Read from ppf file
1349NTuple nt;
1350{
1351PInPersist pis("myfile.ppf");
1352string tag = "NT"; pis.GetObject(nt,tag);
1353}
1354
1355// Fill GeneralData
1356GeneralData mGdata(nt.NEntry());
1357for(int i=0; i<nt.NEntry(); i++)
1358 mGdata.AddData1(xnt[1],xnt[2],xnt[3]); // Fill x, y and error on y
1359mGData.PrintStatus();
1360
1361// Function for fitting : y = f(x) + noise
1362Gauss1DPol mFunction; // gaussian + constant
1363
1364// Prepare for fit
1365GeneralFit mFit(&mFunction); // create a fitter for the choosen function
1366mFit.SetData(&mGData); // connect data to the fitter
1367
1368// Set and initialise the parameters (that's non-linear fitting!)
1369// (num par, name, guess start, step, [limits min and max])
1370mFit.SetParam(0,"high",90.,1..);
1371mFit.SetParam(1,"xcenter",0.05,0.01);
1372mFit.SetParam(2,"sigma",sig,0.05,0.01,10.);
1373 // Give limits to avoid division by zero
1374mFit.SetParam(3,"constant",0.,1.);
1375
1376// Fit and print result
1377int rcfit = mFit.Fit();
1378mFit.PrintFit();
1379if(rcfit>0) {)
1380 cout<<"Reduce_Chisquare = "<<mFit.GetChi2Red()
1381 <<" nstep="<<mFit.GetNStep()<<" rc="<<rcfit<<endl;
1382} else {
1383 cout<<"Fit_Error, rc = "<<rcfit<<" nstep="<<mFit.GetNStep()<<endl;
1384 mFit.PrintFitErr(rcfit);
1385}
1386
1387// Get the result for further use
1388TVector<r_8> ParResult = mFit.GetParm();
1389cout<<ParResult;
1390\end{verbatim}
1391
1392Much more usefull possibilities and detailed informations might be found
1393in the HTML pages of the Sophya manual.
1394
1395\subsection{Polynomial}
1396\index{Polynomial} \index{Poly} \index{Poly2}
1397Polynomials of 1 or 2 variables are supported ({\tt Poly} and {\tt Poly2}).
1398Various operations are supported~:
1399\begin{itemize}
1400\item elementary operations between polynomials $(+,-,*,/) $
1401\item setting or getting coefficients
1402\item computing the value of the polynomial for a given value
1403 of the variable(s),
1404\item derivating
1405\item computing roots (degre 1 or 2)
1406\item fitting the polynomial to vectors of data.
1407\end{itemize}
1408Here is an example of polynomial fitting~:
1409\begin{verbatim}
1410#include "poly.h"
1411// ...
1412Poly pol(2);
1413pol[0] = 100.; pol[1] = 0.; pol[2] = 0.01; // Setting coefficients
1414TVector<r_8> x(100);
1415TVector<r_8> y(100);
1416TVector<r_8> ey(100);
1417for(int i=0;i<100;i++) {
1418 x(i) = i;
1419 ey(i) = 10.;
1420 y(i) = pol((double) i) + ey(i)*NorRand();
1421 ey(i) *= ey(i)
1422}
1423
1424TVector<r_8> errcoef;
1425Poly polfit;
1426polfit.Fit(x,y,ey,2,errcoef);
1427
1428cout<<"Fit Result"<<polfit<<endl;
1429cout<<"Errors :"<<errcoef;
1430\end{verbatim}
1431
1432Similar operations can be done on polynomials with 2 variables.
1433
1434\subsection{Integration, Differential equations}
1435\index{Integration}
1436The NTools module provide also simple classes for numerical integration
1437of functions and differential equations.
1438\begin{figure}[hbt]
1439\dclsbb{Integrator}{GLInteg}
1440\dclsb{TrpzInteg}
1441\end{figure}
1442
1443\index{GLInteg} \index{TrpzInteg}
1444{\bf GLInteg} implements the integration through Gauss-Legendre method
1445and {\bf TrpzInteg} implements trapeze integration. For {\bf TrpzInteg},
1446number of steps specify the number of trapeze, and integration step,
1447their width.
1448The sample code below illustrates the use of TrpzInteg class:
1449\begin{verbatim}
1450#include "integ.h"
1451// ......................................................
1452// Function to be integrated
1453double myf(double x)
1454{
1455// Simple a x + b x^2 (a=2 b=3)
1456return (x*(2.+3.*x));
1457}
1458// ......................................................
1459
1460// Compute Integral(myf, 2., 5.) between xmin=2., xmax=5.
1461TrpzInteg trpz(myf, 2., 5.);
1462// We specify an integration step
1463trpz.DX(0.01);
1464// The integral can be computed as trpz.Value()
1465double myf_integral = trpz.Value();
1466// We could have used the cast operator :
1467cout << "Integral[myf, 2., 5.]= " << (double)trpz << endl;
1468// Limits can be specified through ValueBetween() method
1469cout << "Integral[myf, 0., 4.]= " << trpz.ValueBetween(0.,4.) << endl;
1470\end{verbatim}
1471
1472\subsection{Fourier transform (FFT)}
1473\index{FFT} \index{FFTPackServer}
1474An abstract interface for performing FFT operations is defined by the
1475{\bf FFTServerInterface} class. The {\bf FFTPackSever} class implements
1476one dimensional FFT, on real and complex data. FFTPackServer uses an
1477adapted and extended version of FFTPack (available from netlib),
1478translated in C, and can operate on single and double precision
1479({\tt float, double}) data.
1480
1481The sample code below illustrates the use of FFTServers:
1482\begin{verbatim}
1483#include "fftpserver.h"
1484 // ...
1485TVector<r_8> in(32);
1486TVector< complex<r_8> > out;
1487in = RandomSequence();
1488FFTPackServer ffts;
1489ffts.setNormalize(true); // To have normalized transforms
1490cout << " FFTServer info string= " << ffts.getInfo() << endl;
1491cout << "in= " << in << endl;
1492cout << " Calling ffts.FFTForward(in, out) : " << endl;
1493ffts.FFTForward(in, out);
1494cout << "out= " << out << endl;
1495\end{verbatim}
1496
1497% \newpage
1498\section{Module SUtils}
1499Some utility classes and C/C++ string manipulation functions are gathered
1500in {\bf SUtils} module.
1501\subsection{Using DataCards}
1502\index{DataCards}
1503The {\bf DataCards} class can be used to read parameters from a file.
1504Each line in the file starting with \@ defines a set of values
1505associated with a keyword. In the example below, we read the
1506parameters corresponding with the keyword {\tt SIZE} from the
1507file {\tt ex.d}. We suppose that {\tt ex.d} contains the line: \\
1508{\tt @SIZE 400 250} \\
1509\begin{verbatim}
1510#include "datacards.h"
1511// ...
1512// Initialising DataCards object dc from file ex.d
1513DataCards dc( "ex.d" );
1514// Getting the first and second parameters for keyword size
1515// We define a default value 100
1516int size_x = dc.IParam("SIZE", 0, 100);
1517int size_y = dc.IParam("SIZE", 1, 100);
1518cout << " size_x= " << size_x << " size_y= " << size_y << endl;
1519\end{verbatim}
1520
1521\section{Module SysTools}
1522The {\bf SysTools} module contains classes implementing interface to some
1523OS specific services, such as thread creation and management, dynamic loading and
1524resource usage information. For example, yhe class {\bf Periodic} provides the
1525necessary services needed to implement the execution of a periodic action.
1526
1527\subsection{Resource usage (CPU, memory \ldots) }
1528 The class {\bf ResourceUsage} \index{ResourceUsage}
1529and {\bf Timer} \index{Timer} provides access to information
1530about various resource usage (memory, CPU, ...).
1531The class {\bf Timer} \index{time (CPU, elapsed)} and c-functions
1532{\tt InitTim() , PrtTim(const char * Comm) } can be used to print
1533the amount of CPU and elapsed time in programs.
1534
1535The following sample code illustrates the use of {\bf ResourceUsage} :
1536\begin{verbatim}
1537 // How to check resource usage for a given part of the program
1538 ResourceUsage res;
1539 // --- Part of the program to be checked : Start
1540 // ...
1541 res.Update();
1542 cout << " Memory size increase (KB):" << res.getDeltaMemorySize() << endl;
1543 cout << " Resource usage info : \n" << res << endl;
1544\end{verbatim}
1545
1546\subsection{Thread management classes}
1547\index{ZThread} \index{ZMutex}
1548A basic interface to POSIX threads is also provided
1549through the \index{threads} {\bf ZThread}, {\bf ZMutex} and {\bf ZSync}
1550classes. The best way to use thread management classes is by inheriting
1551from {\bf ZThread} and redefining the {\tt run() } method.
1552It is also possible to use the default {\tt run() } implementation and associate
1553a function to perform the action, as in the example below :
1554\begin{verbatim}
1555 // The functions to perform computing
1556 void fun1(void *arg) { }
1557 void fun2(void *arg) { }
1558 // ...
1559 ZThread zt1;
1560 zt1.setAction(fun1, arg[1]);
1561 ZThread zt2;
1562 zt2.setAction(fun2, arg[1]);
1563 cout << " Starting threads ... " << endl;
1564 zt1.start();
1565 zt2.start();
1566 cout << " Waiting for threads to end ... " << endl;
1567 zt1.join();
1568 zt2.join();
1569\end{verbatim}
1570The classes {\bf ZMutex} \index{mutex} and {\bf ZSync} can be used
1571to perform synchronisation and signaling between threads.
1572
1573\subsection{Dynamic linker and C++ compiler classes}
1574\index{PDynLinkMgr}
1575The class {\bf PDynLinkMgr} can be used for managing shared libraries
1576at run time. The example below shows the run time linking of a function:\\
1577{\tt extern "C" { void myfunc(); } } \\
1578\begin{verbatim}
1579#include "pdlmgr.h"
1580// ...
1581string soname = "mylib.so";
1582string funcname = "myfunc";
1583PDynLinkMgr dyl(soname);
1584DlFunction f = dyl.GetFunction(funcname);
1585if (f != NULL) {
1586// Calling the function
1587 f();
1588}
1589\end{verbatim}
1590
1591\index{CxxCompilerLinker}
1592The {\bf CxxCompilerLinker} class provides the services to compile C++ code and building
1593shared libraries, using the same compiler and options which have
1594been used to create the SOPHYA shared library.
1595The sample program below illustrates using this class to build
1596the shared library (myfunc.so) from the source file myfunc.cc :
1597\begin{verbatim}
1598#include "cxxcmplnk.h"
1599// ...
1600string flnm = "myfunc.cc";
1601string oname, soname;
1602int rc;
1603CxxCompilerLinker cxx;
1604// The Compile method provides a default object file name
1605rc = cxx.Compile(flnm, oname);
1606if (rc != 0 ) { // Error when compiling ... }
1607// The BuildSO method provides a default shared object file name
1608rc = cxx.BuildSO(oname, soname);
1609if (rc != 0 ) { // Error when creating shared object ... }
1610\end{verbatim}
1611
1612\subsection{Command interpreter}
1613The class {\bf Commander} can be used in interactive programs to provide
1614c-shell like command interpreter and scripting capabilties.
1615Arithmetic expression evaluation is implemented through the {\bf CExpressionEvaluator}
1616and {\bf RPNExpressionEvaluator} classes.
1617The command language provides variable manipulation through the usual
1618{\tt \$varname} vector variable and arithmetic expression extensions, as well
1619as the control and test blocs.
1620\begin{verbatim}
1621#include "commander.h"
1622...
1623Commander cmd;
1624char* ss[3] = {"foreach f ( AA bbb CCCC ddddd )", "echo $f" , "end"};
1625for(int k=0; k<3; k++) {
1626 string line = ss[k];
1627 cmd.Interpret(line);
1628}
1629\end{verbatim}
1630
1631\newpage
1632\section{Module SkyMap}
1633\begin{figure}[hbt]
1634\dclsbb{AnyDataObj}{PixelMap}
1635\dclsccc{PixelMap}{Sphericalmap}{SphereHEALPix}
1636\dclsc{SphereThetaPhi}
1637\dclsc{SphereECP}
1638\dclsb{LocalMap}
1639\caption{partial class diagram for spherical map classes in Sophya}
1640\end{figure}
1641The {\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.
1642\subsection {Spherical maps}
1643SkyMap module provides three kinds of complete ($4 \pi$) spherical maps according to the
1644pixelization scheme.
1645SphereHEALPix represents spheres pixelized following the HEALPIix algorithm (E. Hivon, K. Gorski)
1646\footnote{see the HEALPix Homepage: http://www.eso.org/kgorski/healpix/ }
1647, 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) :
1648\index{\tcls{SphereHEALPix}}
1649
1650\begin{verbatim}
1651#include "spherehealpix.h"
1652// ...
1653SphereHEALPix<double> sph(8);
1654for (int k=0; k< sph.NbPixels(); k++) sph(k) = (double)(10*k);
1655\end{verbatim}
1656
1657SphereThetaPhi is used in a similar way with an argument representing number of slices in theta (Euler angle) for an hemisphere.
1658\index{\tcls{SphereThetaPhi}}
1659The SphereECP class correspond to the cylindrical projection and can be used for representing
1660partial or full spherical maps. However, it has the disadvantage of having non uniform pixel
1661size.
1662\index{\tcls{SphereECP}}
1663
1664\subsection {Local maps}
1665\index{\tcls{LocalMap}}
1666A 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).
1667
1668Internally, 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(...))
1669
1670The 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).
1671\begin{verbatim}
1672#include "localmap.h"
1673//..............
1674 LocalMap<r_4> locmap(4,5);
1675 for (int k=0; k<locmap.NbPixels();k++) locmap(k)=10.*k;
1676 locmap.SetOrigin();
1677 locmap.SetSize(30.,30.);
1678\end{verbatim}
1679
1680\subsection{Writing, viewing \dots }
1681
1682All these objects have been design to be written to or read from a persistant file.
1683The following example shows how to write the previously created objects
1684into such a file~:
1685\begin{verbatim}
1686//-- Writing
1687
1688#include "fiospherehealpix.h"
1689//................
1690
1691char *fileout = "myfile.ppf";
1692POutPersist outppf(fileout);
1693FIO_SphereHEALPix<r_8> outsph(sph);
1694outsph.Write(outppf);
1695FIO_LocalMap<r_8> outloc(locmap);
1696outloc.Write(outppf);
1697// It is also possible to use the << operator
1698POutPersist os("sph.ppf");
1699os << outsph;
1700os << outloc;
1701\end{verbatim}
1702
1703Sophya graphical tools (spiapp) can automatically display and operate
1704all these objects.
1705
1706\newpage
1707\section{Samba and SkyT}
1708\subsection{Samba}
1709\index{Spherical Harmonics}
1710\index{SphericalTransformServer}
1711The 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...
1712\begin{verbatim}
1713#include "skymap.h"
1714#include "samba.h"
1715....................
1716
1717// Generate input spectra a + b* l + c * gaussienne(l, 50, 20)
1718int lmax = 92;
1719Vector clin(lmax);
1720for(int l=0; l<lmax; l++) {
1721 double xx = (l-50.)/10.;
1722 clin(l) = 1.e-2 -1.e-4*l + 0.1*exp(-xx*xx);
1723}
1724
1725// Compute map from spectra
1726SphericalTransformServer<r_8> ylmserver;
1727int m = 128; // HealPix pixelisation parameter
1728SphereHEALPix<r_8> map(m);
1729ylmserver.GenerateFromCl(map, m, clin, 0.);
1730// Compute power spectrum from map
1731Vector clout = ylmserver.DecomposeToCl(map, lmax, 0.);
1732\end{verbatim}
1733
1734\subsection{Module SkyT}
1735\index{RadSpectra} \index{SpectralResponse}
1736The SkyT module is composed of two types of classes:
1737\begin{itemize}
1738\item{} one which corresponds to an emission spectrum of
1739radiation, which is called RadSpectra
1740\item{} one which corresponds to the spectral response
1741of a given detector (i.e. corresponding to a detector
1742filter in a given frequency domain), which is called
1743SpectralResponse.
1744\end{itemize}
1745\begin{figure}[hbt]
1746\dclsbb{RadSpectra}{RadSpectraVec}
1747\dclsb{BlackBody}
1748\dclsccc{AnyDataObj}{SpectralResponse}{SpecRespVec}
1749\dclsc{GaussianFilter}
1750\caption{partial class for SkyT module}
1751\end{figure}
1752
1753\begin{verbatim}
1754#include "skyt.h"
1755// ....
1756// Compute the flux from a blackbody at 2.73 K through a square filter
1757BlackBody myBB(2.73);
1758// We define a square filter from 100 - 200 GHz
1759SquareFilter mySF(100,200);
1760// Compute the filtered integrated flux :
1761double flux = myBB.filteredIntegratedFlux(mySF);
1762\end{verbatim}
1763
1764A more detailed description of SkyT module can be found in:
1765{\it The SkyMixer (SkyT and PMixer modules) - Sophya Note No 2. }
1766available also from Sophya Web site.
1767
1768\newpage
1769\section{Module FitsIOServer}
1770This module provides classes for handling file input-output in FITS
1771\footnote{http://heasarc.gsfc.nasa.gov/docs/software/fitsio/fitsio.html}
1772format using the cfitsio library. Its
1773design is similar to the SOPHYA persistence (see Module BaseTools).
1774Delegate classes or handlers perform the actual read/write from/to fits files.
1775\par
1776Compared to the SOPHYA native persistence (PPF format),
1777FITS format has the advantage of being used extensively, and handled
1778by a many different software tools. It is a de facto standard in
1779astronomy and astrophysics.
1780However, FITS lacks some of the features present in SOPHYA PPF, and although
1781many SOPHYA objects can be saved in FITS files, FITS persistence has
1782some limitations. For example, FITS does not currently handle complex arrays.
1783\subsection{FITS streams}
1784\index{FITS} \index{FitsInOutFile}
1785%%
1786The class {\bf FitsInOutFile} can be seen as a wrapper class for the cfitsio library functions.
1787This class has been introduced in 2005 (V=1.9), when the module has been
1788extensively changed. In order to keep backward compatibility, the old fits wrapper
1789classes ({\bf FitsFile, FitsInFile, FitsOutFile}) has been changed to inherit from
1790 {\bf FitsInOutFile}. The use of class FitsFile and specific services of these old classes
1791 should be avoided, but FitsInFile, FitsOutFile can be safely considered a specialisation
1792 of FitsInOutFile for read/input and write/output operations respectively.
1793 Most c-fitsio errors are converted to an exception: {\bf FitsIOException}.
1794 \par
1795 File names are passed to cfitsio library. It is thus possible to use cfitsio file name conventions,
1796 such as {\bf ! } as the first character, for overwriting existing files (when creating files).
1797 The diagram below shows the class hierarchy for cfitsio wrapper classes.
1798\begin{figure}[hbt]
1799\dclsa{FitsInOutFile}
1800\dclscc{FitsFile}{FitsInFile}
1801\dclscc{FitsFile}{FitsOutFile}
1802\end{figure}
1803%%%%
1804\subsection{FITS handlers}
1805\index{FitsManager}
1806Handlers classes inheriting from {\bf FitsHandlerInterface} perform write/read operations
1807for AnyDataObj objects to/from FitsInOutFile streams. The {\bf FitsManager} class provides
1808top level services for object read/write in FITS files.
1809\par In most cases,
1810\hspace{5mm} {\tt FitsInOutFile\& $<<$ } \, and \, {\tt FitsInOutFile\& $>>$ } \hspace{5mm}
1811operators can be used to write and read objects.
1812\par
1813The two main types of fits data structures, images and tables
1814{\tt (IMAGE\_HDU , BINARY\_TBL , ASCII\_TBL)} are handled by the generic handlers: \\
1815{\bf \tcls{FitsArrayHandler}} and {\bf FitsHandler$<$BaseDataTable$>$}.
1816\par
1817A number of more specific handlers are also available, in particular for NTuple,
1818\tcls{SphereHealPix} and \tcls{SphereThetaPhi}. \\[2mm]
1819{\bf Warning:} Some handlers were written with the old FitsIOServer classes and
1820inherit from the intermediate class {\bf FitsIOHandler}. They
1821have been adapted to the new scheme, but may not perform correctly if not used
1822through the {\bf FitsManager} class. \\[2mm]
1823%%%
1824The examples below illustrates the usage of FitsIOServer classes. They can be compiled
1825and executed using runcxx, without the {\tt include} lines: \\[1mm]
1826\hspace*{5mm} {\tt csh> runcxx -import SkyMap -import FitsIOServer -inc fiosinit.h }
1827%%%
1828\begin{enumerate}
1829\item Saving an array and a HealPix map to a Fits file
1830\begin{verbatim}
1831#include "fitsioserver.h"
1832#include "fiosinit.h"
1833// ....
1834{
1835// Make sure FitsIOServer module is initialised :
1836FitsIOServerInit();
1837// Create and open a fits file named myfile.fits
1838FitsInOutFile fos("myfile.fits", FitsInOutFile ::Fits_Create);
1839// Create and save a 15x11 matrix of integers
1840TMatrix<int_4> mxi(15, 11);
1841mxi = RegularSequence(10.,10.);
1842fos << mxi;
1843// Save a HEALPix spherical map using FitsManager services
1844SphereHEALPix<r_8> sph(16);
1845sph = 48.3;
1846FitsManager::Write(fos, sph);
1847// --- The << operator could have been used instead : fos << sph;
1848}
1849\end{verbatim}
1850%%%%
1851%%%%
1852\item Reading objects and the header from the previously created fits file:
1853\begin{verbatim}
1854{
1855FitsIOServerInit(); // Initialisation
1856// ---- Open the fits file named myfile.fits
1857FitsInFile fis("myfile.fits");
1858//---- print file information on cout
1859cout << fis << endl;
1860//--- Read in the array
1861TArray<int_4> arr;
1862fis >> arr;
1863arr.Show();
1864//--- Position on second HDU
1865fis.MoveAbsToHDU(2);
1866//--- read and display header information
1867DVList hdu2;
1868fis.GetHeaderRecords(hdu2, true, true);
1869cout << hdu2;
1870//--- read in the HEALPix map
1871SphereHEALPix<r_8> sph;
1872FitsManager::Read(fis, sph);
1873// --- The >> operator could have been used instead : fis >> sph;
1874sph.Show();
1875}
1876\end{verbatim}
1877%%%%%%%
1878%%%
1879\item DataTable objects can be read from and written to FITS files as ASCII or
1880binary tables. The example belo show reading the DataTable created in the example
1881in section \ref{datatables} from a PPF file and saving it to a fits file.
1882\begin{verbatim}
1883#include "swfitsdtable.h"
1884// ....
1885{
1886FitsIOServerInit(); // FitsIOServer Initialisation
1887//--- Reading in DataTable object from PPF file
1888PInPersist pin("dtable.ppf");
1889DataTable dt;
1890pin >> dt;
1891dt.Show();
1892//--- Saving table to FITS
1893FitsInOutFile fos("!dtable.fits", FitsInOutFile ::Fits_Create);
1894fos << dt;
1895}
1896\end{verbatim}
1897%%%%
1898\end{enumerate}
1899%%%
1900A partial class diagram of FITS persistence handling classes is shown below. The
1901class {\tt FitsIOhandler} conforms to the old FitsIOServer module design and
1902should not be used anymore.
1903\begin{figure}[hbt]
1904\dclsbb{FitsHandlerInterface}{FitsArrayHandler$<$T$>$}
1905\dclsb{\tcls{FitsHandler}}
1906\dclscc{FitsIOhandler}{FITS\_NTuple}
1907\dclsc{FITS\_SphereHEALPix}
1908% \dclsb{FITS\_LocalMap}
1909\end{figure}
1910
1911\subsection{SwFitsDataTable and other classes}
1912\index{SwFitsDataTable}
1913The {\bf SwFitsDataTable} class implements the BaseDataTable interface
1914using a FITS file as swap space. Compared to SwPPFDataTable, they can be
1915used in R/W mode (reading from the table, when it is being created / filled).
1916They can be used in a way similar to DataTable and SwPPFDataTable.
1917When creating the table, a {\tt FitsInOutFile } stream, opened for writing has
1918to be passed to the creator. No further operation is needed.
1919\begin{verbatim}
1920// ....
1921FitsInOutFile so("!myswtable.fits", FitsInOutFile::Fits_Create);
1922SwFitsDataTable dt(so, 16);
1923// define table columns
1924dt.AddIntegerColumn("X0_i");
1925dt.AddFloatColumn("X1_f");
1926// ... Fill the table
1927r_8 x[5];
1928for(int i=0; i<63; i++) {
1929 x[0] = (i%9)-4.; x[1] = (i/9)-3.;
1930 dt.AddLine(x);
1931}
1932\end{verbatim}
1933The class {\bf FitsBTNtuIntf } provide an alternative tool to read FITS tables.
1934{\bf FitsABTColRd} , {\bf FitsABTWriter } and {\bf FitsImg2DWriter } can also
1935be used to manipulate FITS files.
1936\par
1937The {\bf scanfits} program can be used to check FITS files and analyse their
1938content (See \ref{scanfits}).
1939
1940%%%%
1941\newpage
1942\section{LinAlg and IFFTW modules}
1943An interface to use LAPACK library (available from {\tt http://www.netlib.org})
1944is implemented by the {\bf LapackServer} class, in module LinAlg.
1945\index{LapackServer}.
1946The sample code below shows how to use SVD (Singular Value Decomposition)
1947through LapackServer:
1948\begin{verbatim}
1949#include "intflapack.h"
1950// ...
1951// Use FortranMemoryMapping as default
1952BaseArray::SetDefaultMemoryMapping(BaseArray::FortranMemoryMapping);
1953// Create an fill the arrays A and its copy AA
1954int n = 20;
1955Matrix A(n , n), AA;
1956A = RandomSequence(RandomSequence::Gaussian, 0., 4.);
1957AA = A; // AA is a copy of A
1958// Compute the SVD decomposition
1959Vector S; // Vector of singular values
1960Matrix U, VT;
1961LapackServer<r_8> lpks;
1962lpks.SVD(AA, S, U, VT);
1963// We create a diagonal matrix using S
1964Matrix SM(n, n);
1965for(int k=0; k<n; k++) SM(k,k) = S(k);
1966// Check the result : A = U*SM*VT
1967Matrix diff = U*(SM*VT) - A;
1968double min, max;
1969diff.MinMax(min, max);
1970cout << " Min/Max difference Matrix (?=0) , Min= " << min
1971 << " Max= " << max << endl;
1972\end{verbatim}
1973
1974\index{FFTWServer}
1975The {\bf FFTWServer} class (in module FFTW) implements FFTServerInterface class
1976methods, for one dimensional and multi-dimensional Fourier
1977transforms on double precision data using the FFTW package
1978(available from {\tt http://www.fftw.org}).
1979
1980\newpage
1981\section{Building and installing Sophya}
1982\subsection{supported platforms}
1983Presently, the Sophya library has been tested with the following
1984compiler/platform pairs:
1985
1986\begin{center}
1987\begin{tabular}{|l|l|}
1988\hline
1989OS & compiler \\
1990\hline
1991Linux (RH) & g++ (3.2) \\
1992Linux (SCL) & icc (8.1) (Intel compiler) \\
1993MacOSX/Darwin 10.3 & g++ 3.3 \\
1994HP/Compaq/DEC Tru64 ( OSF1) & cxx (6.1 , 6.3) \\
1995SGI IRIX64 & CC (7.3) \\
1996IBM AIX & xlC (7.x) \\
1997\hline
1998\end{tabular}
1999\end{center}
2000
2001Some of the modules in the Sophya package uses external libraries. The
2002{\bf FitsIOServer} is the example of such a module, where the {\tt libcfitsio.a}
2003is used.
2004\par
2005The object files from a given Sophya module are grouped in an archive library
2006with the module's name ({\tt libmodulename.a}). All Sophya modules
2007 are grouped in a single shared library ({\tt libsophya.so}), while the
2008modules with reference to external libraries are grouped in
2009({\tt libextsophya.so}). The {\bf PI} and {\bf PIext} modules are
2010grouped in ({\tt libPI.so}).
2011Alternatively, it is possible to group all modules in a single shared
2012library {\tt libAsophyaextPI.so}.
2013\par
2014Each library module has a {\tt Makefile} which compiles the source files
2015and build the correspond static (archive) library using the compilation
2016rules and flags defined in {\tt \$SOPHYABASE/include/sophyamake.inc}.
2017Each program module has a {\tt Makefile} which compiles and link the
2018corresponding programs using the compilation rules and libraries
2019defined in {\$SOPHYABASE/include/sophyamake.inc}.
2020The top level Makefile in BuildMgr/ compiles each library modules
2021and builds shared libraries.
2022\par
2023The series of Makefiles use the link to {\tt sophyamake.inc} in BuildMgr.
2024There are also the {\tt smakefile} series which uses the explicit path, using
2025{\tt \$SOPHYABASE} environment variable.
2026
2027\subsection{Installation}
2028\label{build}
2029The build procedure has two main steps:
2030\begin{enumerate}
2031\item The configure step (BuildMgr/configure) setup the directory structure and
2032the necessary configuration file. Refer to section \ref{directories} for
2033the description of SOPHYA directory tree and files.
2034\item The make step compiles the different sources files, create the library and optionaly
2035builds all or some of the associated executables.
2036\end{enumerate}
2037
2038{\tt BuildMgr/configure } is a c-shell script with a number of arguments:
2039\begin{verbatim}
2040csh> ./configure -h
2041configure [-sbase SOPHYABASE] [-scxx SOPHYACXX] [-incln]
2042 [-minc mymake.inc]
2043 [-extp dir1 -extp dir2 ...] [-extip dir1 -extip dir2 ... ]
2044 [-extlp dir1 -extlp dir2 ... ]
2045 [-noextlib -noext fits -noext fftw -noext lapack ]
2046 [-noext astro -noext minuit]
2047 [-usefftw2 -uselapack2] [-singleslb]
2048\end{verbatim}
2049\begin{itemize}
2050\item[] -sbase : define SOPHYA installation base directory. \$SOPHYABASE is used
2051if not specified.
2052\item[] -scxx : selects the C++ compiler. \$SOPHYACXX s used
2053if not specified.
2054\item[] -incln : creates symbolic link for include files, instead of copying them.
2055\item[] -minc : give an explicit name for the file used to generate
2056\$SOPHYABASE/include/sophyamake.inc.
2057\item[] -extp : Adds the specied path to the search path of the external libraries
2058include files and archive library.
2059\item[] -extip : Adds the specied path to the search path of the external libraries
2060include files.
2061\item[] -extp : Adds the specied path to the search path of the external libraries
2062archive (libxxx.a).
2063\item[] -noextlib : Disable compiling of modules referencing external libraries.
2064\item[] -noext : Disable compiling of the specified module (with reference to external
2065library.
2066\item[] -usefftw2: FFTW V2 is being used (default FFTW V3) - A compilation flag
2067will be defined in sspvflags.h
2068\item[] -uselapack2: Lapack V2 is being used (defaulr V3) - A compilation flag
2069will be defined in sspvflags.h
2070\item[] -singleslb: A single shared library for all SOPHYA, PI and external library interface
2071modules will be build. A compilation flag
2072will be defined in sspvflags.h . `See also target {\tt slballinone} below.
2073\end{itemize}
2074
2075In the example below, we assume that we want to install Sophya from a
2076released (tagged) version in the source directory {\tt \$SRC} in the
2077{\tt /usr/local/Sophya} directory, using {\tt g++}. We assume that
2078the external libraries can be found in {\tt /usr/local/ExtLibs/}.
2079We disable the compilation of the MinuitAdapt and XAstrPack packages.
2080
2081\vspace*{3mm}
2082\begin{verbatim}
2083# Create the top level directory
2084csh> mkdir /usr/local/Sophya/
2085csh> cd $SRC/BuildMgr
2086# Step 1.a : Run the configuration script
2087csh> ./configure -sbase /usr/local/Sophya -scxx g++ -extp /usr/local/ExtLibs/ \
2088-noext astro -noext minuit
2089# Step 1.b : Check the generated file $SOPHYABASE/include/sophyamake.inc
2090csh> ls -lt *.inc
2091csh> more sophyamake.inc
2092\end{verbatim}
2093If necessary, edit the generated file {\tt sophyamake.inc } in order to modify
2094compilation flags, library list. The file is rather short and self documented.
2095\begin{verbatim}
2096# Step 2.a: Compile the modules without external library reference
2097csh> make libs
2098# Step 2.b: Compile the modules WITH external library reference (optional)
2099csh> make extlibs
2100# Step 2.c: Build libsophya.so
2101csh> make slb
2102# Step 2.d: Build libextsophya.so (optional)
2103csh> make slbext
2104# Step 2.e: Compile the PI and PIext modules (optional)
2105csh> make PI
2106# Step 2.f: Build the corresponding shared library libPI.so (optional)
2107csh> make slbpi
2108\end{verbatim}
2109
2110To compile all modules and build the shared libraries, it is possible
2111to perform the steps 2.a to 2.f using the targets {\tt all} and {\tt slball}
2112defined in the Makefile
2113\begin{verbatim}
2114# Step 2.a ... 2.f
2115csh> make all slball
2116\end{verbatim}
2117
2118It is also possible to group all modules in a single shared library using
2119the target {\tt slballinone}.
2120\begin{verbatim}
2121# Step 2.a ... 2.f
2122csh> make all slballinone
2123\end{verbatim}
2124
2125At this step, all libraries should have been made. Programs using
2126Sophya libraries can now be built:
2127\begin{verbatim}
2128# To compile some of the test programs
2129csh> make basetests
2130# To compile runcxx , scanppf , scanfits
2131csh> make prgutil
2132# To build (s)piapp (libPI.so is needed)
2133csh> make piapp
2134\end{verbatim}
2135
2136If no further modification or update of source files is foreseen, it is possible
2137to remove all .o files:
2138\begin{verbatim}
2139# To clean $SOPHYABASE/obj directory :
2140csh> make cleanobj
2141\end{verbatim}
2142
2143
2144\subsection{Notes}
2145\begin{itemize}
2146\item[{\bf Makefile}] List of top level Makefile build targets
2147\begin{verbatim}
2148> libs extlibs PI = all
2149> slb slbext slbpi = slball (OR = slballinone)
2150> clean cleanobj
2151> tests prgutil prgmap progpi = prgall
2152> basetests piapp (ou progpi) pmixer
2153\end{verbatim}
2154\item[{\bf MacOS X}] A high performance mathematic and signal processing
2155library, including LAPACK and BLAS is packaged in Darwin/MacOS X (10.3, 10.4) : \\
2156\hspace*{5mm} {\bf -framework Accelerate}
2157\item[{\bf Tru64/OSF}] An optimised math library with LAPACK and BLAS might
2158optionaly be installed {\bf (-lcxlm) }. On our system, this libray contained Lapack V2.
2159So we used the LAPACK, as compiled from the public sources, and linked with
2160the Tru64 native BLAS.
2161\item[{\bf IRIX64}] We used the math library with LAPACK V2 and BLAS
2162from SGI : {\bf -lcomplib.sgimath}
2163\item[{\bf AIX}] There seem to be a problem on AIX when several shared
2164libraries are used. We have been able to run SOPHYA programs either
2165using static libraries, or a single shared library (libAsophyaextPI.so)
2166if extlibs and PI are needed, in addition to stand alone SOPHYA modules.
2167It has not been possible to link SOPHYA with fortran libraries
2168\item[{\bf Mgr}] This module contains makefiles and build scripts
2169that were used in SOPHYA up to version 1.7 (2004) : OBSOLETE.
2170\end{itemize}
2171
2172\subsection{Files and scripts in BuildMgr/ }
2173\begin{itemize}
2174\item[] {\bf Makefile:} Top level Makefile for building SOPHYA.
2175{\tt smakefile} is similar to Makefile, except that it uses
2176the smakefiles in each module.
2177\item[] {\bf mkmflib:} c-shell script for creation of library module
2178Makefile / smakefile.
2179{\tt ./mkmflib -sbase /tmp/sbase SUtils }
2180\item[] {\b mkmfprog:}
2181{\tt c-shell script for creation of programs module
2182Makefile / smakefile \\
2183{\tt ./mkmfprog -sbase /tmp/sbase ProgPI }
2184\item[] {\bf domkmf:} c-shell script - calls mkmflib for all modules \\
2185{\tt ./domkmf -sbase /tmp/sbase}
2186\item[] Configuration files for different compilers and OS
2187(\tt Linux\_g++\_make.inc , OSF1\_cxx\_make.inc \ldots }.
2188These files are used to generate {\tt sophyamake.inc}
2189\end{itemize}
2190
2191
2192
2193\newpage
2194\appendix
2195\section{SOPHYA Exceptions}
2196\index{Exception classes} \index{PThrowable} \index{PError} \index{PException}
2197SOPHYA library defines a set of exceptions which are used
2198for signalling error conditions. The figure below shows a partial
2199class diagram for exception classes in SOPHYA.
2200\begin{figure}[hbt]
2201\dclsbb{PThrowable}{PError}
2202\dclscc{PError}{AllocationError}
2203\dclscc{PError}{NullPtrError}
2204\dclscc{PError}{ForbiddenError}
2205\dclscc{PError}{AssertionFailedError}
2206\dclsbb{PThrowable}{PException}
2207\dclscc{PException}{IOExc}
2208\dclscc{PException}{SzMismatchError}
2209\dclscc{PException}{RangeCheckError}
2210\dclscc{PException}{ParmError}
2211\dclscc{PException}{TypeMismatchExc}
2212\dclscc{PException}{MathExc}
2213\dclscc{PException}{CaughtSignalExc}
2214\caption{partial class diagram for exception handling in Sophya}
2215\end{figure}
2216
2217For simple programs, it is a good practice to handle
2218the exceptions at least at high level, in the {\tt main()} function.
2219The example below shows the exception handling and the usage
2220of Sophya persistence.
2221
2222\input{ex1.inc}
2223
2224
2225\newpage
2226\addcontentsline{toc}{section}{Index}
2227\printindex
2228\end{document}
2229
Note: See TracBrowser for help on using the repository browser.