Array.h

Go to the documentation of this file.

//# Array.h: A templated N-D Array class with zero origin
//# Copyright (C) 1993,1994,1995,1996,1997,1998,1999,2000,2001,2002,2003
//# Associated Universities, Inc. Washington DC, USA,
//#
//# This library is free software; you can redistribute it and/or modify it
//# under the terms of the GNU Library General Public License as published by
//# the Free Software Foundation; either version 2 of the License, or (at your
//# option) any later version.
//#
//# This library is distributed in the hope that it will be useful, but WITHOUT
//# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
//# FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Library General Public
//# License for more details.
//#
//# You should have received a copy of the GNU Library General Public License
//# along with this library; if not, write to the Free Software Foundation,
//# Inc., 675 Massachusetts Ave, Cambridge, MA 02139, USA.
//#
//# Correspondence concerning AIPS++ should be addressed as follows:
//#        Internet email: aips2-request@nrao.edu.
//#        Postal address: AIPS++ Project Office
//#                        National Radio Astronomy Observatory
//#                        520 Edgemont Road
//#                        Charlottesville, VA 22903-2475 USA
//#
//# $Id: Array.h 20739 2009-09-29 01:15:15Z Malte.Marquarding $

#ifndef CASA_ARRAY_H
#define CASA_ARRAY_H

//# Includes
#include <casa/Arrays/ArrayBase.h>
#include <casa/Containers/Block.h>
#include <casa/Utilities/CountedPtr.h>
#include <casa/Arrays/MaskLogiArrFwd.h>
#include <casa/Arrays/IPosition.h>
#include <casa/ostream.h>
#include <iterator>
#if defined(WHATEVER_VECTOR_FORWARD_DEC)
WHATEVER_VECTOR_FORWARD_DEC;
#else
#include <casa/stdvector.h>
#endif

namespace casa { //#Begin casa namespace
//# Forward Declarations
class AipsIO;
class Slice;
class Slicer;
template<class T> class Matrix;
template<class T> class ArrayIterator;
template<class T> class MaskedArray;
template<class Domain, class Range> class Functional;
//template <class T, class U> class vector; 


// <summary> A templated N-D Array class with zero origin </summary>

// Array<T> is a templated, N-dimensional, Array class. The origin is zero,
// but by default indices are zero-based. This Array class is the
// base class for specialized Vector<T>, Matrix<T>, and Cube<T> classes.
//
// Indexing into the array, and positions in general, are given with IPosition
// (essentially a vector of integers) objects. That is, an N-dimensional 
// array requires a length-N IPosition to define a position within the array.
// Unlike C, indexing is done with (), not []. Also, the storage order
// is the same as in FORTRAN, i.e. memory varies most rapidly with the first
// index.
// <srcblock>
//                                     // axisLengths = [1,2,3,4,5]
// IPosition axisLengths(5, 1, 2, 3, 4, 5); 
// Array<Int> ai(axisLengths);         // ai is a 5 dimensional array of
//                                     // integers; indices are 0-based
//                                     // => ai.nelements() == 120
// Array<Int> ai2(axisLengths);        // The first element is at index 0
// IPosition zero(5); zero = 0;        // [0,0,0,0,0]
// //...
// </srcblock>
// Indexing into an N-dimensional array is relatively expensive. Normally
// you will index into a Vector, Matrix, or Cube. These may be obtained from
// an N-dimensional array by creating a reference, or by using an 
// ArrayIterator. The "shape" of the array is an IPosition which gives the
// length of each axis.
//
// An Array may be standalone, or it may refer to another array, or to
// part of another array (by refer we mean that if you change a pixel in 
// the current array, a pixel in the referred to array also changes, i.e.
// they share underlying storage).
// <note role=warning>
//        One way one array can reference another is through the copy 
//        constructor. While this might be what you want, you should
//        probably use the reference() member function to make it explicit.
//        The copy constructor is used when arguments are passed by value;
//        normally functions should not pass Arrays by value, rather they
//        should pass a reference or a const reference. On the positive
//        side, returning an array from a function is efficient since no
//        copying need be done.
// </note>
//
// Aside from the explicit reference() member function, a user will
// most commonly encounter an array which references another array
// when he takes an array slice (or section). A slice is a sub-region of
// an array (which might also have a stride: every nth row, every mth column,
// ...).
// <srcblock>
// IPosition lengths(3,10,20,30);
// Array<Int> ai(lengths);         // A 10x20x30 cube
// Cube<Int> ci;
// //...
// ci.reference(ai1);              // ci and ai now reference the same
//                                 // storage
// ci(0,0,0) = 123;                // Can use Cube indexing
// ci.xyPlane(2) = 0;              //     and other member functions
// IPosition zero(3,0,0,0);
// assert(ai(zero) == 123);        // True because ai, ci are references
// //...
// Array<Int> subArray;
// IPosition blc(3,0,0,0), trc(3,5,5,5);
// subArray.reference(ai(blc, trc));
// subArray = 10;                  // All of subArray, which is the
//                                 // subcube from 0,0,0 to 5,5,5 in
//                                 // ai, has the value 10.
// </srcblock>
// While the last example has an array slice referenced explicitly by another
// array variable, normally the user will often only use the slice as
// a temporary in an expresion, for example:
// <srcblock>
//   Array<Complex> array;
//   IPosition blc, trc, offset;
//   //...
//   // Copy from one region of the array into another
//   array(blc, trc) = array(blc+offset, trc+offset);
// </srcblock>
//
// The Array classes are intended to operate on relatively large
// amounts of data. While they haven't been extensively tuned yet,
// they are relatively efficient in terms of speed. Presently they
// are not space efficient -- the overhead is about 15 words. While
// this will be improved (probably to about 1/2 that), these array
// classes are not appropriate for very large numbers of very small
// arrays. The Block<T> class may be what you want in this circumstance.
//
// Element by element mathematical and logical operations are available
// for arrays (defined in aips/ArrayMath.h and aips/ArrayLogical.h).
// Because arithmetic and logical functions are split out, it is possible
// to create an Array<T> (and hence Vector<T> etc) for any type T that has
// a default constructor, assignment operator, and copy constructor. In
// particular, Array<String> works.
//
// If compiled with the preprocessor symbol AIPS_DEBUG symbol, array
// consistency ("invariants") will be checked in most member
// functions, and indexing will be range-checked. This should not be
// defined for production runs.
//
// A tutorial for the ArrayClasses is available in the "AIPS++ Programming
// Manual."
//
// <note role=tip>
// Most of the data members and functions which are "protected" should
// likely become "private".
// </note>
//
// <todo asof="1999/12/30">
//   <li> Integrate into the Lattice hierarchy
//   <li> Factor out the common functions (shape etc) into a type-independent
//        base class.
// </todo>

template<class T> class Array : public ArrayBase
{
public:

    // Result has dimensionality of zero, and  nelements is zero.
    Array();

    // Create an array of the given shape, i.e. after construction
    // array.ndim() == shape.nelements() and array.shape() == shape.
    // The origin of the Array is zero.
    explicit Array(const IPosition &shape);

    // Create an array of the given shape and initialize it with the
    // initial value.
    Array(const IPosition &shape, const T &initialValue);

    // After construction, this and other reference the same storage.
    Array(const Array<T> &other);

    // Create an Array of a given shape from a pointer.
    Array(const IPosition &shape, T *storage, StorageInitPolicy policy = COPY);
    // Create an Array of a given shape from a pointer. Because the pointer
    // is const, a copy is always made.
    Array(const IPosition &shape, const T *storage);

    // Frees up storage only if this array was the last reference to it.
    virtual ~Array();

    // Assign the other array to this array.
    // If the shapes mismatch, this array is resized.
    virtual void assign (const Array<T>& other);

    // Set every element of the array to "value." Also could use the
    // assignment operator which assigns an array from a scalar.
    void set(const T &value);

    // Apply the function to every element of the array. This modifies
    // the array in place.
    // <group>
    // This version takes a function which takes a T and returns a T.
    void apply(T (*function)(T));
    // This version takes a function which takes a const T reference and
    // returns a T.
    void apply(T (*function)(const T &));
    // This version applies a functional.
    void apply(const Functional<T,T> &function);
    // </group>

    // After invocation, this array and other reference the same storage. That
    // is, modifying an element through one will show up in the other. The
    // arrays appear to be identical; they have the same shape.
    // <br>Please note that this function makes it possible to reference a
    // const Array, thus effectively it makes a const Array non-const.
    // Although this may seem undesirable at first sight, it is necessary to
    // be able to make references to temporary Array objects, in particular to 
    // Array slices. Otherwise one first needs to use the copy constructor.
    //# The const has been introduced on 2005-Mar-31 because of the hassle
    //# involved in calling the copy ctor before reference.
    virtual void reference(const Array<T> &other);

    // Copy the values in other to this. If the array on the left hand
    // side has no elements, then it is resized to be the same size as
    // as the array on the right hand side. Otherwise, the arrays must
    // conform (same shapes).
    // <srcblock>
    // IPosition shape(2,10,10);     // some shape
    // Array<Double> ad(shape);
    // //...
    // Array<Double> ad2;            // N.B. ad2.nelements() == 0
    // ad2 = ad;                     // ad2 resizes, then elements
    //                               //     are copied.
    // shape = 20;
    // Array<Double> ad3(shape);
    // ad3 = ad;                     // Error: arrays do not conform
    // </srcblock>
    // Note that the assign function can be used to assign a
    // non-conforming array.
    virtual Array<T>