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

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

Suite preparation documentation SOPHYA pour V2.0 , Reza 17/07/2006

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