PLearn 0.1
Tensor.h
Go to the documentation of this file.
00001 // -*- C++ -*-
00002 
00003 // TTensor.h: Definition of a tensor
00004 // Copyright (c) 2002 by Julien Keable and Pascal Vincent
00005 
00006 // Redistribution and use in source and binary forms, with or without
00007 // modification, are permitted provided that the following conditions are met:
00008 // 
00009 //  1. Redistributions of source code must retain the above copyright
00010 //     notice, this list of conditions and the following disclaimer.
00011 // 
00012 //  2. Redistributions in binary form must reproduce the above copyright
00013 //     notice, this list of conditions and the following disclaimer in the
00014 //     documentation and/or other materials provided with the distribution.
00015 // 
00016 //  3. The name of the authors may not be used to endorse or promote
00017 //     products derived from this software without specific prior written
00018 //     permission.
00019 // 
00020 // THIS SOFTWARE IS PROVIDED BY THE AUTHORS ``AS IS'' AND ANY EXPRESS OR
00021 // IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
00022 // OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN
00023 // NO EVENT SHALL THE AUTHORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
00024 // SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED
00025 // TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00026 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00027 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00028 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00029 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00030 // 
00031 // This file is part of the PLearn library. For more information on the PLearn
00032 // library, go to the PLearn Web site at www.plearn.org
00033 
00034 #ifndef TENSOR_INC
00035 #define TENSOR_INC
00036 
00037 #include <iostream.h>
00038 #include <plearn/base/Storage.h>
00039 #include <plearn/base/TinyVector.h>
00040 #include <vector>
00041 
00042 using namespace std;
00043 
00044 namespace PLearn {
00045 
00046 template <class T>
00047 class TTensor;
00048 
00049 typedef TinyVector<int, 7> IVec;
00050 
00051 /*
00052   class: TTensorElementIterator 
00053   Iterator that iterates over all elements of a TTensor
00054 */
00055 
00056 template<class T>
00057 class TTensorElementIterator
00058 {
00059 private:
00060     // tensor from which we iterate
00061     TTensor<T>* tensor;
00062     IVec position_;
00063     // the element 'i' in the following vector is the memory leap from the last elem. of dimension i to the next element 
00064     // computed on creation to accelerate computations in operator++
00065     IVec stride_minus_width;
00066     T* ptr; // current element pointer
00067     // is the iterator past the last tensor element?
00068     bool end_met;
00069 
00070 public:
00071 
00072     inline TTensorElementIterator()
00073         :tensor(NULL),ptr(NULL),end_met(false)
00074     {
00075     }
00076 
00077     inline TTensorElementIterator(TTensor<T>* tensor_, const IVec& pos)
00078         :tensor(tensor_),position_(pos),end_met(false)
00079     {
00080         ptr = tensor->data()+tensor->linearIndex(pos);
00081         stride_minus_width.resize(tensor->ndims());
00082         for(int i=0;i<(signed)tensor->ndims()-1;i++)
00083             stride_minus_width[i] = tensor->stride_[i+1] - tensor->stride_[i] * tensor->width_[i];
00084     }
00085 
00086     inline TTensorElementIterator<T>& operator++();
00087 
00088     inline T* operator->() const
00089     { return ptr; }
00090   
00091     inline T& operator*() const
00092     { return *ptr; }
00093 
00094     inline bool operator==(const TTensorElementIterator& other)
00095     { return ptr==other.ptr; }
00096 
00097     inline bool operator!=(const TTensorElementIterator& other)
00098     { return ptr!=other.ptr; }
00099     
00100     inline IVec position() const {return position_;}
00101 
00102     inline bool end() const {return end_met;}
00103 };
00104 
00105 /*
00106   class TTensorSubTensorIterator :
00107   this iterator iterates on subTensor of a parent Tensor.
00108   For exemple, you might have a 3d tensor (x,y,z) and iterate over all (x,y) planes.
00109   To build such an object, call TTensor::getSubTensorIterator(IVec) with a vector
00110   containing non-zero values for the dimensions that the subTensors spread on
00111 */
00112 
00113 template<class T>
00114 class TTensorSubTensorIterator
00115 {
00116 private:
00117     // tensor from which we iterate
00118     TTensor<T>* tensor;
00119     IVec position_;
00120     IVec dimensions_;
00121     // is the iterator past the last element?
00122     bool end_met;
00123 
00124 public:
00125 
00126     inline TTensorSubTensorIterator(TTensor<T>* tensor_, const IVec& dim)
00127         :tensor(tensor_),dimensions_(dim),end_met(false)
00128     {
00129         position_=IVec((unsigned)tensor->ndims(),0);
00130         for(int i=0;i<tensor->ndims();i++)
00131             if(dimensions_[i]!=0)
00132                 position_[i]=-1;
00133     }
00134 
00135     inline TTensorSubTensorIterator<T>()
00136         :tensor(NULL),end_met(false){}
00137 
00138     inline IVec position() const {return position_;}
00139 
00140     inline TTensorSubTensorIterator<T>& operator++();
00141 
00142     inline TTensor<T> operator*() const
00143     {
00144         IVec from(tensor->ndims()),length(tensor->ndims());
00145         for(int i=0;i<tensor->ndims();i++)
00146         {
00147             if(position_[i]==-1)
00148             {
00149                 from[i]=0;
00150                 length[i]=tensor->width_[i];
00151             }
00152             else
00153             {
00154                 from[i]=position_[i];
00155                 length[i]=1;
00156             }
00157         }
00158         return tensor->subTensor(from,length);
00159     }
00160 
00161     inline bool end() const {return end_met;}
00162 };
00163 
00164 
00165 
00166 /*
00167   class Tensor: 
00168   A tensor is a n-dimension matrix. The number of dimensions is determined
00169   at run-time upon object creation. The tensor, as of today, cannot be resized 
00170   (that needs to be implemented, the default constructor is currently pointless)
00171   You can extract a subTensor by specifying two vectors : one for the start 
00172   indices and one for the widths of the subTensor on each dimensions.
00173   Dimensions with width == 0 are not 'dimensions' anymore. E.g: a(x,y) plane 
00174   extracted from a 3d tensor (x,y,z) has only 2 dimensions.
00175 */
00176  
00177 template <class T>
00178 class TTensor
00179 {
00180     friend class TTensorElementIterator<T>;
00181     friend class TTensorSubTensorIterator<T>;
00182 protected:
00183     int offset_; 
00184     IVec stride_; 
00185     IVec width_; 
00186     PP< Storage<T> > storage; 
00188 public:
00189     typedef TTensorElementIterator<T> iterator; // iterator over elements
00190     typedef TTensorSubTensorIterator<T> subTensorIterator; // iterator over elements
00191     
00192     TTensor<T>()
00193         :offset_(0), stride_(0), width_(0)
00194     {}
00195 
00196     TTensor<T>(IVec len)
00197         :offset_(0), stride_(0), width_(len)
00198     { 
00199         stride_.resize(len.size());
00200         for(int i=0;i<ndims();i++)
00201         {
00202             stride_[i]=1;
00203             for(int j=i-1;j>=0;j--)
00204                 stride_[i]*=len[j];
00205         }
00206         storage = new Storage<T>(totalElements());
00207     }
00208 
00209     TTensor<T>(const IVec& len, const T& init_value)
00210         :offset_(0), stride_(0), width_(len)
00211     { 
00212         stride_.resize(len.size());
00213         for(int i=0;i<ndims();i++)
00214         {
00215             stride_[i]=1;
00216             for(int j=i-1;j>=0;j--)
00217                 stride_[i]*=len[j];
00218         }
00219         storage = new Storage<T>(totalElements());
00220         fill(init_value);
00221     }
00222     
00224     int resize(const IVec& len)
00225     {
00226         width_=len;
00227         stride_.resize(len.size());
00228         for(int i=0;i<ndims();i++)
00229         {
00230             stride_[i]=1;
00231             for(int j=i-1;j>=0;j--)
00232                 stride_[i]*=len[j];
00233         }
00234         storage = new Storage<T>(totalElements());
00235         offset_=0;
00236     }
00237     
00238     int ndims() const { return width_.size();}
00239 
00240     IVec width() const {return width_;}
00241     IVec sizes() const {return width_;}
00242     int size(int k) const { return width_[k]; }
00243 
00244     int totalElements() const 
00245     {
00246         if(ndims()==0)
00247             return 0;
00248         int te=1;
00249         for(int i=0;i<ndims();i++)
00250             te*=width_[i];
00251         return te;
00252     }
00253 
00254     int linearIndex(const IVec& pos) const
00255     {
00256         int idx=0;
00257         for(int i=0;i<ndims();i++)
00258             idx+=pos[i]*stride_[i];
00259         return idx;
00260     }
00261     
00262     int linearIndex(const vector<int>& pos) const
00263     {
00264         int idx=0;
00265         for(int i=0;i<ndims();i++)
00266             idx+=pos[i]*stride_[i];
00267         return idx;
00268     }
00269     
00270       
00271     // extract a sub tensor from this one. 'from' is a vector of start indices in each dimension, len is the length of each new dimension
00272     // Any dimension of length 0 will not be considered as a dimension anymore (we'll keep a slice at the particular position given by from).
00273     // i.e. : a subTensor extracted from a 3d tensor (x,y,z) with the len of the z dimension set to 0 will be a 2d tensor
00274     TTensor<T> subTensor(const IVec& from, const IVec& len)
00275     {
00276         TTensor<T> subt = *this;
00277         subt.width_ = len;
00278         subt.offset_ = linearIndex(from);
00279         
00280         IVec idx;
00281         for(int i=0;i<ndims();i++)
00282         {
00283             if(from[i]<0 || from[i]+len[i]>width_[i] || len[i]<0)
00284                 PLERROR("TTensor::subTensor : at index %i : from, len, width = %i %i %i",i,from[i],len[i],width_[i]);
00285             if(len[i]>0) // skip the 0 dimensions
00286                 idx.push_back(i);      
00287         }
00288         // if idx is empty, it is because all lengths are 0: we have a scalar
00289         if(idx.empty())
00290             idx.push_back(0);
00291         subt.selectDimensions(idx);
00292         return subt;
00293     }
00294 
00295 
00296 
00297     // NOTE: (Pascal) This call has been deprecated because the throw_useless_dimension behaviour is unsatisfactory and a potential for nasty bugs. 
00298     // I've rewritten another subTensor method with a clearer semantic: a len of 1 means keep that dimension (as a single element with that particular value)
00299     // and a len of 0 means throw away that dimension.
00300     // extract a sub tensor from this one. 'from' is a vector of start indices in each dimension, len is the length of each new dimension
00301     // ** NOTE :
00302     // If 'throw_useless_dimensions' == true, any dimension of length 1 will not be considered as a dimension anymore
00303     // i.e. : a subTensor extracted from a 3d matrix (x,y,z) with a fixed z value will be a 2d matrix rather than a 3d matrix 
00304     // with a width of one in the 'z' dimension
00305     TTensor<T> DEPRECATEDsubTensor(const IVec& from, const IVec& len, bool throw_useless_dimensions=true)
00306     {
00307         TTensor<T> subt = *this;
00308         subt.width_ = len;
00309         subt.offset_ = linearIndex(from);
00310         
00311         IVec idx;
00312         for(int i=0;i<ndims();i++)
00313         {
00314             if(from[i]<0 || from[i]+len[i]>width_[i] || len[i]<0)
00315                 PLERROR("TTensor::subTensor : at index %i : from, len, width = %i %i %i",i,from[i],len[i],width_[i]);
00316             if(len[i]>1 || !throw_useless_dimensions)
00317                 idx.push_back(i);      
00318         }
00319         // if idx is empty, it is because all lengths are 1 since a single element was selected,
00320         // so we need to create a single dimension.
00321         if(idx.empty())
00322             idx.push_back(0);
00323         subt.selectDimensions(idx);
00324         return subt;
00325     }
00326 
00327     TTensor<T> operator[](int i)
00328     {
00329         IVec from;
00330         IVec len;
00331         from.push_back(i);
00332         len.push_back(0);
00333         for(int k=1; k<ndims(); k++)
00334         {
00335             from.push_back(0);
00336             len.push_back(size(k));
00337         }
00338         return subTensor(from,len);
00339     }
00340 
00341     void selectDimensions(const IVec& dim)
00342     {
00343         IVec newwidth,newstride;
00344         for(int i=0;i<(signed)dim.size();i++)
00345         {
00346             newwidth.push_back(width_[dim[i]]);
00347             newstride.push_back(stride_[dim[i]]);
00348         }
00349         stride_=newstride;
00350         width_=newwidth;
00351     }
00352    
00353     T& operator()(const IVec& pos) const {return (*storage)[linearIndex(pos)];}
00354     T& operator()(const vector<int>& pos) const {return (*storage)[linearIndex(pos)];}
00355 
00356     void fill(const T& val)
00357     {
00358         int te=totalElements();
00359         for(int i=0;i<te;i++)
00360             (*storage)[i]=val;
00361     }
00362     
00363     inline T* data() const
00364     {
00365 #ifdef BOUNDCHECK
00366         if(!storage)
00367             PLERROR("IN TTensor::data()\nAttempted to get a pointer to the data of an empty matrix");
00368 #endif
00369         return storage->data+offset_; 
00370     }
00371 
00372     IVec lastElementPos() const
00373     {
00374         IVec v;
00375         v.resize(ndims());
00376         for(int i=0;i<ndims();i++)
00377             v[i]=(width_[i]-1);
00378         return v;
00379     }
00380     
00382     inline iterator begin()
00383     { 
00384         return iterator(this, IVec((unsigned)ndims(),0));
00385     }
00386     
00387     // note that Tensor iterators know if they're past the end of data by themselves
00388     // this method is only implemented to provide a standard interface. Check for iterator.end()
00389     inline iterator end()
00390     { 
00391         return iterator(this,lastElementPos());
00392     }
00393     
00395     inline subTensorIterator getSubTensorIterator(const IVec& v)
00396     { 
00397         return subTensorIterator(this, v);
00398     }
00399 
00400 };
00401 
00402 typedef TTensor<real> Tensor;
00403 
00404 template<class T>
00405 inline TTensorElementIterator<T>& TTensorElementIterator<T>::operator++()
00406 { 
00407     ptr+=tensor->stride_[0];
00408     position_[0]++;
00409     // if we hit the end of the first dimension
00410     // we need to compute a few things
00411     if(position_[0]==tensor->width_[0])
00412     {
00413         bool found=false;
00414         int idx=0;
00415         // find the lowest index of the position vector we can increment
00416         // a.k.a get the next element of the tensor
00417         while(!found && idx<(signed int)(tensor->ndims()-1))
00418         {
00419             idx++;
00420             found = position_[idx]<tensor->width_[idx]-1;
00421         }
00422         if(found)
00423         {
00424             position_[idx]++;
00425             for(int i=0;i<idx;i++)
00426             {
00427                 ptr+=stride_minus_width[i];
00428                 position_[i]=0;
00429             }      
00430         }
00431         else 
00432             // all elements of the tensor have been visited           
00433             end_met=true;
00434     }
00435     return *this;
00436 }
00437     
00438 
00439 template<class T>
00440 inline TTensorSubTensorIterator<T>& TTensorSubTensorIterator<T>::operator++()
00441 {
00442     bool found=false;
00443     int idx=-1;
00444     // find the lowest index of the position vector we can increment
00445     while(!found && idx<(signed int)(tensor->ndims()-1))
00446     {
00447         idx++;
00448         found = position_[idx]!=-1 && position_[idx]<tensor->width_[idx]-1;
00449     }
00450     if(found)
00451     {
00452         position_[idx]++;
00453         for(int i=0;i<idx;i++)
00454             if(position_[i]!=-1)
00455                 position_[i]=0;
00456     }
00457     else 
00458         // all subTensors of the tensor have been visited           
00459         end_met=true;
00460     
00461     return *this;
00462 }
00463 
00464  
00465 };
00466 
00467 #endif
00468 
00469 
00470 /*
00471   Local Variables:
00472   mode:c++
00473   c-basic-offset:4
00474   c-file-style:"stroustrup"
00475   c-file-offsets:((innamespace . 0)(inline-open . 0))
00476   indent-tabs-mode:nil
00477   fill-column:79
00478   End:
00479 */
00480 // vim: filetype=cpp:expandtab:shiftwidth=4:tabstop=8:softtabstop=4:encoding=utf-8:textwidth=79 :
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines