subview_cube_meat.hpp

Go to the documentation of this file.
00001 // Copyright (C) 2009 NICTA
00002 // 
00003 // Authors:
00004 // - Conrad Sanderson (conradsand at ieee dot org)
00005 // 
00006 // This file is part of the Armadillo C++ library.
00007 // It is provided without any warranty of fitness
00008 // for any purpose. You can redistribute this file
00009 // and/or modify it under the terms of the GNU
00010 // Lesser General Public License (LGPL) as published
00011 // by the Free Software Foundation, either version 3
00012 // of the License or (at your option) any later version.
00013 // (see http://www.opensource.org/licenses for more info)
00014 
00015 
00016 //! \addtogroup subview_cube
00017 //! @{
00018 
00019 
00020 template<typename eT>
00021 inline
00022 subview_cube<eT>::~subview_cube()
00023   {
00024   arma_extra_debug_sigprint();
00025   }
00026 
00027 
00028 
00029 template<typename eT>
00030 arma_inline
00031 subview_cube<eT>::subview_cube
00032   (
00033   const Cube<eT>& in_m,
00034   const u32       in_row1,
00035   const u32       in_col1,
00036   const u32       in_slice1,
00037   const u32       in_row2,
00038   const u32       in_col2,
00039   const u32       in_slice2
00040   )
00041   : m           (in_m)
00042   , m_ptr       (0)
00043   , aux_row1    (in_row1)
00044   , aux_col1    (in_col1)
00045   , aux_slice1  (in_slice1)
00046   , aux_row2    (in_row2)
00047   , aux_col2    (in_col2)
00048   , aux_slice2  (in_slice2)
00049   , n_rows      (1 + in_row2 - in_row1)
00050   , n_cols      (1 + in_col2 - in_col1)
00051   , n_elem_slice(n_rows * n_cols)
00052   , n_slices    (1 + in_slice2 - in_slice1)
00053   , n_elem      (n_elem_slice * n_slices)
00054   {
00055   arma_extra_debug_sigprint();
00056   }
00057 
00058 
00059 
00060 template<typename eT>
00061 arma_inline
00062 subview_cube<eT>::subview_cube
00063   (
00064         Cube<eT>& in_m,
00065   const u32       in_row1,
00066   const u32       in_col1,
00067   const u32       in_slice1,
00068   const u32       in_row2,
00069   const u32       in_col2,
00070   const u32       in_slice2
00071   )
00072   : m           (in_m)
00073   , m_ptr       (&in_m)
00074   , aux_row1    (in_row1)
00075   , aux_col1    (in_col1)
00076   , aux_slice1  (in_slice1)
00077   , aux_row2    (in_row2)
00078   , aux_col2    (in_col2)
00079   , aux_slice2  (in_slice2)
00080   , n_rows      (1 + in_row2 - in_row1)
00081   , n_cols      (1 + in_col2 - in_col1)
00082   , n_elem_slice(n_rows * n_cols)
00083   , n_slices    (1 + in_slice2 - in_slice1)
00084   , n_elem      (n_elem_slice * n_slices)
00085   {
00086   arma_extra_debug_sigprint();
00087   }
00088 
00089 
00090 
00091 template<typename eT>
00092 inline
00093 void
00094 subview_cube<eT>::operator+= (const eT val)
00095   {
00096   arma_extra_debug_sigprint();
00097   
00098   for(u32 slice = 0; slice < n_slices; ++slice)
00099     {
00100     for(u32 col = 0; col < n_cols; ++col)
00101       {
00102       eT* coldata = slice_colptr(slice,col);
00103       
00104       for(u32 row = 0; row < n_rows; ++row)
00105         {
00106         coldata[row] += val;
00107         }
00108       
00109       }
00110     }
00111 
00112   }
00113 
00114 
00115 
00116 template<typename eT>
00117 inline
00118 void
00119 subview_cube<eT>::operator-= (const eT val)
00120   {
00121   arma_extra_debug_sigprint();
00122   
00123   for(u32 slice = 0; slice < n_slices; ++slice)
00124     {
00125     for(u32 col = 0; col<n_cols; ++col)
00126       {
00127       eT* coldata = slice_colptr(slice,col);
00128       
00129       for(u32 row = 0; row<n_rows; ++row)
00130         {
00131         coldata[row] -= val;
00132         }
00133       
00134       }
00135     }
00136   }
00137 
00138 
00139 
00140 template<typename eT>
00141 inline
00142 void
00143 subview_cube<eT>::operator*= (const eT val)
00144   {
00145   arma_extra_debug_sigprint();
00146   
00147   for(u32 slice = 0; slice < n_slices; ++slice)
00148     {
00149     for(u32 col = 0; col<n_cols; ++col)
00150       {
00151       eT* coldata = slice_colptr(slice,col);
00152       
00153       for(u32 row = 0; row<n_rows; ++row)
00154         {
00155         coldata[row] *= val;
00156         }
00157       
00158       }
00159     }  
00160   }
00161 
00162 
00163 
00164 template<typename eT>
00165 inline
00166 void
00167 subview_cube<eT>::operator/= (const eT val)
00168   {
00169   arma_extra_debug_sigprint();
00170   
00171   for(u32 slice = 0; slice < n_slices; ++slice)
00172     {
00173     for(u32 col = 0; col<n_cols; ++col)
00174       {
00175       eT* coldata = slice_colptr(slice,col);
00176       
00177       for(u32 row = 0; row<n_rows; ++row)
00178         {
00179         coldata[row] /= val;
00180         }
00181       
00182       }
00183     }  
00184   }
00185 
00186 
00187 
00188 template<typename eT>
00189 template<typename T1>
00190 inline
00191 void
00192 subview_cube<eT>::operator= (const BaseCube<eT,T1>& in)
00193   {
00194   arma_extra_debug_sigprint();
00195   
00196   const unwrap_cube<T1> tmp(in.get_ref());
00197   
00198   const Cube<eT>&         x = tmp.M;
00199         subview_cube<eT>& t = *this;
00200   
00201   arma_debug_assert_same_size(t, x, "copy into subcube");
00202   
00203   
00204   for(u32 slice = 0; slice < t.n_slices; ++slice)
00205     {
00206     for(u32 col = 0; col < t.n_cols; ++col)
00207       {
00208             eT* t_coldata = t.slice_colptr(slice,col);
00209       const eT* x_coldata = x.slice_colptr(slice,col);
00210       
00211       for(u32 row = 0; row < t.n_rows; ++row)
00212         {
00213         t_coldata[row] = x_coldata[row];
00214         }
00215         
00216       }
00217     }
00218   }
00219 
00220 
00221 
00222 template<typename eT>
00223 template<typename T1>
00224 inline
00225 void
00226 subview_cube<eT>::operator+= (const BaseCube<eT,T1>& in)
00227   {
00228   arma_extra_debug_sigprint();
00229   
00230   const unwrap_cube<T1> tmp(in.get_ref());
00231   
00232   const Cube<eT>&         x = tmp.M;
00233         subview_cube<eT>& t = *this;
00234   
00235   arma_debug_assert_same_size(t, x, "cube addition");
00236   
00237   for(u32 slice = 0; slice < t.n_slices; ++slice)
00238     {
00239     for(u32 col = 0; col < t.n_cols; ++col)
00240       {
00241             eT* t_coldata = t.slice_colptr(slice,col);
00242       const eT* x_coldata = x.slice_colptr(slice,col);
00243       
00244       for(u32 row = 0; row < t.n_rows; ++row)
00245         {
00246         t_coldata[row] += x_coldata[row];
00247         }
00248       }
00249     }
00250   
00251   }
00252 
00253 
00254 
00255 template<typename eT>
00256 template<typename T1>
00257 inline
00258 void
00259 subview_cube<eT>::operator-= (const BaseCube<eT,T1>& in)
00260   {
00261   arma_extra_debug_sigprint();
00262   
00263   const unwrap_cube<T1> tmp(in.get_ref());
00264   
00265   const Cube<eT>&         x = tmp.M;
00266         subview_cube<eT>& t = *this;
00267   
00268   arma_debug_assert_same_size(t, x, "cube subtraction");
00269   
00270   
00271   for(u32 slice = 0; slice < t.n_slices; ++slice)
00272     {
00273     for(u32 col = 0; col < t.n_cols; ++col)
00274       {
00275             eT* t_coldata = t.slice_colptr(slice,col);
00276       const eT* x_coldata = x.slice_colptr(slice,col);
00277       
00278       for(u32 row = 0; row < t.n_rows; ++row)
00279         {
00280         t_coldata[row] -= x_coldata[row];
00281         }
00282       }
00283     }  
00284   }
00285 
00286 
00287 
00288 template<typename eT>
00289 template<typename T1>
00290 inline
00291 void
00292 subview_cube<eT>::operator%= (const BaseCube<eT,T1>& in)
00293   {
00294   arma_extra_debug_sigprint();
00295   
00296   const unwrap_cube<T1> tmp(in.get_ref());
00297   
00298   const Cube<eT>&         x = tmp.M;
00299         subview_cube<eT>& t = *this;
00300   
00301   arma_debug_assert_same_size(t, x, "cube schur product");
00302   
00303   for(u32 slice = 0; slice < t.n_slices; ++slice)
00304     {
00305     for(u32 col = 0; col<t.n_cols; ++col)
00306       {
00307             eT* t_coldata = t.slice_colptr(slice,col);
00308       const eT* x_coldata = x.slice_colptr(slice,col);
00309       
00310       for(u32 row = 0; row<t.n_rows; ++row)
00311         {
00312         t_coldata[row] *= x_coldata[row];
00313         }
00314       }
00315     }
00316   
00317   }
00318 
00319 
00320 
00321 template<typename eT>
00322 template<typename T1>
00323 inline
00324 void
00325 subview_cube<eT>::operator/= (const BaseCube<eT,T1>& in)
00326   {
00327   arma_extra_debug_sigprint();
00328   
00329   const unwrap_cube<T1> tmp(in.get_ref());
00330   
00331   const Cube<eT>&         x = tmp.M;
00332         subview_cube<eT>& t = *this;
00333   
00334   arma_debug_assert_same_size(t, x, "element-wise cube division");
00335   
00336   
00337   for(u32 slice = 0; slice < t.n_slices; ++slice)
00338     {
00339     for(u32 col = 0; col<t.n_cols; ++col)
00340       {
00341             eT* t_coldata = t.slice_colptr(slice,col);
00342       const eT* x_coldata = x.slice_colptr(slice,col);
00343       
00344       for(u32 row = 0; row<t.n_rows; ++row)
00345         {
00346         t_coldata[row] /= x_coldata[row];
00347         }
00348       }
00349     }
00350   
00351   }
00352 
00353 
00354 
00355 //! x.subcube(...) = y.subcube(...)
00356 template<typename eT>
00357 inline
00358 void
00359 subview_cube<eT>::operator= (const subview_cube<eT>& x_in)
00360   {
00361   arma_extra_debug_sigprint();
00362   
00363   const bool overlap = check_overlap(x_in);
00364   
00365         Cube<eT>*         tmp_cube         = overlap ? new Cube<eT>(x_in.m) : 0;
00366   const subview_cube<eT>* tmp_subview_cube = overlap ? new subview_cube<eT>(*tmp_cube, x_in.aux_row1, x_in.aux_col1, x_in.aux_slice1, x_in.aux_row2, x_in.aux_col2, x_in.aux_slice2) : 0;
00367   const subview_cube<eT>&                x = overlap ? (*tmp_subview_cube) : x_in;
00368   
00369   subview_cube<eT>& t = *this;
00370   
00371   arma_debug_assert_same_size(t, x, "copy into subcube");
00372   
00373   
00374   for(u32 slice = 0; slice < t.n_slices; ++slice)
00375     {
00376     for(u32 col = 0; col < t.n_cols; ++col)
00377       {
00378             eT* t_coldata = t.slice_colptr(slice,col);
00379       const eT* x_coldata = x.slice_colptr(slice,col);
00380       
00381       for(u32 row = 0; row < t.n_rows; ++row)
00382         {
00383         t_coldata[row] = x_coldata[row];
00384         }
00385       }
00386     }
00387     
00388   if(overlap)
00389     {
00390     delete tmp_subview_cube;
00391     delete tmp_cube;
00392     }
00393   
00394   }
00395 
00396 
00397 
00398 template<typename eT>
00399 inline
00400 void
00401 subview_cube<eT>::operator+= (const subview_cube<eT>& x_in)
00402   {
00403   arma_extra_debug_sigprint();
00404   
00405   const bool overlap = check_overlap(x_in);
00406   
00407         Cube<eT>*         tmp_cube         = overlap ? new Cube<eT>(x_in.m) : 0;
00408   const subview_cube<eT>* tmp_subview_cube = overlap ? new subview_cube<eT>(*tmp_cube, x_in.aux_row1, x_in.aux_col1, x_in.aux_slice1, x_in.aux_row2, x_in.aux_col2, x_in.aux_slice2) : 0;
00409   const subview_cube<eT>&                x = overlap ? (*tmp_subview_cube) : x_in;
00410   
00411   subview_cube<eT>& t = *this;
00412   
00413   arma_debug_assert_same_size(t, x, "cube addition");
00414   
00415   
00416   for(u32 slice = 0; slice < t.n_slices; ++slice)
00417     {
00418     for(u32 col = 0; col < t.n_cols; ++col)
00419       {
00420             eT* t_coldata = t.slice_colptr(slice,col);
00421       const eT* x_coldata = x.slice_colptr(slice,col);
00422       
00423       for(u32 row = 0; row < t.n_rows; ++row)
00424         {
00425         t_coldata[row] += x_coldata[row];
00426         }
00427       }
00428     }
00429     
00430   if(overlap)
00431     {
00432     delete tmp_subview_cube;
00433     delete tmp_cube;
00434     }
00435   
00436   }
00437 
00438 
00439 
00440 template<typename eT>
00441 inline
00442 void
00443 subview_cube<eT>::operator-= (const subview_cube<eT>& x_in)
00444   {
00445   arma_extra_debug_sigprint();
00446   
00447   const bool overlap = check_overlap(x_in);
00448   
00449         Cube<eT>*         tmp_cube         = overlap ? new Cube<eT>(x_in.m) : 0;
00450   const subview_cube<eT>* tmp_subview_cube = overlap ? new subview_cube<eT>(*tmp_cube, x_in.aux_row1, x_in.aux_col1, x_in.aux_slice1, x_in.aux_row2, x_in.aux_col2, x_in.aux_slice2) : 0;
00451   const subview_cube<eT>&                x = overlap ? (*tmp_subview_cube) : x_in;
00452   
00453   subview_cube<eT>& t = *this;
00454   
00455   arma_debug_assert_same_size(t, x, "cube subtraction");
00456   
00457   
00458   for(u32 slice = 0; slice < t.n_slices; ++slice)
00459     {
00460     for(u32 col = 0; col < t.n_cols; ++col)
00461       {
00462             eT* t_coldata = t.slice_colptr(slice,col);
00463       const eT* x_coldata = x.slice_colptr(slice,col);
00464       
00465       for(u32 row = 0; row < t.n_rows; ++row)
00466         {
00467         t_coldata[row] -= x_coldata[row];
00468         }
00469       }
00470     }
00471     
00472   if(overlap)
00473     {
00474     delete tmp_subview_cube;
00475     delete tmp_cube;
00476     }
00477     
00478   }
00479 
00480 
00481 
00482 template<typename eT>
00483 inline
00484 void
00485 subview_cube<eT>::operator%= (const subview_cube& x_in)
00486   {
00487   arma_extra_debug_sigprint();
00488   
00489   const bool overlap = check_overlap(x_in);
00490   
00491         Cube<eT>*         tmp_cube         = overlap ? new Cube<eT>(x_in.m) : 0;
00492   const subview_cube<eT>* tmp_subview_cube = overlap ? new subview_cube<eT>(*tmp_cube, x_in.aux_row1, x_in.aux_col1, x_in.aux_slice1, x_in.aux_row2, x_in.aux_col2, x_in.aux_slice2) : 0;
00493   const subview_cube<eT>&                x = overlap ? (*tmp_subview_cube) : x_in;
00494   
00495   subview_cube<eT>& t = *this;
00496   
00497   arma_debug_assert_same_size(t, x, "element-wise cube multiplication");
00498   
00499   
00500   for(u32 slice = 0; slice < t.n_slices; ++slice)
00501     {
00502     for(u32 col = 0; col < t.n_cols; ++col)
00503       {
00504             eT* t_coldata = t.slice_colptr(slice,col);
00505       const eT* x_coldata = x.slice_colptr(slice,col);
00506       
00507       for(u32 row = 0; row < t.n_rows; ++row)
00508         {
00509         t_coldata[row] *= x_coldata[row];
00510         }
00511       }
00512     }
00513     
00514   if(overlap)
00515     {
00516     delete tmp_subview_cube;
00517     delete tmp_cube;
00518     }
00519   
00520   }
00521 
00522 
00523 
00524 template<typename eT>
00525 inline
00526 void
00527 subview_cube<eT>::operator/= (const subview_cube& x_in)
00528   {
00529   arma_extra_debug_sigprint();
00530   
00531   const bool overlap = check_overlap(x_in);
00532   
00533         Cube<eT>*         tmp_cube         = overlap ? new Cube<eT>(x_in.m) : 0;
00534   const subview_cube<eT>* tmp_subview_cube = overlap ? new subview_cube<eT>(*tmp_cube, x_in.aux_row1, x_in.aux_col1, x_in.aux_slice1, x_in.aux_row2, x_in.aux_col2, x_in.aux_slice2) : 0;
00535   const subview_cube<eT>&                x = overlap ? (*tmp_subview_cube) : x_in;
00536   
00537   subview_cube<eT>& t = *this;
00538   
00539   arma_debug_assert_same_size(t, x, "element-wise cube division");
00540   
00541   
00542   for(u32 slice = 0; slice < t.n_slices; ++slice)
00543     {
00544     for(u32 col = 0; col < t.n_cols; ++col)
00545       {
00546             eT* t_coldata = t.slice_colptr(slice,col);
00547       const eT* x_coldata = x.slice_colptr(slice,col);
00548       
00549       for(u32 row = 0; row < t.n_rows; ++row)
00550         {
00551         t_coldata[row] /= x_coldata[row];
00552         }
00553       }
00554     }
00555     
00556   if(overlap)
00557     {
00558     delete tmp_subview_cube;
00559     delete tmp_cube;
00560     }
00561   
00562   }
00563 
00564 
00565 
00566 template<typename eT>
00567 template<typename T1>
00568 inline
00569 void
00570 subview_cube<eT>::operator= (const Base<eT,T1>& in)
00571   {
00572   arma_extra_debug_sigprint();
00573   
00574   const unwrap<T1> tmp(in.get_ref());
00575   
00576   const Mat<eT>&          x = tmp.M;
00577         subview_cube<eT>& t = *this;
00578   
00579   arma_debug_assert_same_size(t, x, "copy into subcube");
00580   
00581   
00582   for(u32 col = 0; col < t.n_cols; ++col)
00583     {
00584           eT* t_coldata = t.slice_colptr(t.aux_slice1, col);
00585     const eT* x_coldata = x.colptr(col);
00586     
00587     for(u32 row = 0; row < t.n_rows; ++row)
00588       {
00589       t_coldata[row] = x_coldata[row];
00590       }
00591       
00592     }
00593   }
00594 
00595 
00596 
00597 template<typename eT>
00598 template<typename T1>
00599 inline
00600 void
00601 subview_cube<eT>::operator+= (const Base<eT,T1>& in)
00602   {
00603   arma_extra_debug_sigprint();
00604   
00605   const unwrap<T1> tmp(in.get_ref());
00606   
00607   const Mat<eT>&          x = tmp.M;
00608         subview_cube<eT>& t = *this;
00609   
00610   arma_debug_assert_same_size(t, x, "cube addition");
00611   
00612   for(u32 col = 0; col < t.n_cols; ++col)
00613     {
00614           eT* t_coldata = t.slice_colptr(t.aux_slice1, col);
00615     const eT* x_coldata = x.colptr(col);
00616     
00617     for(u32 row = 0; row < t.n_rows; ++row)
00618       {
00619       t_coldata[row] += x_coldata[row];
00620       }
00621     }
00622   }
00623 
00624 
00625 
00626 template<typename eT>
00627 template<typename T1>
00628 inline
00629 void
00630 subview_cube<eT>::operator-= (const Base<eT,T1>& in)
00631   {
00632   arma_extra_debug_sigprint();
00633   
00634   const unwrap<T1> tmp(in.get_ref());
00635   
00636   const Mat<eT>&          x = tmp.M;
00637         subview_cube<eT>& t = *this;
00638   
00639   arma_debug_assert_same_size(t, x, "cube subtraction");
00640   
00641   for(u32 col = 0; col < t.n_cols; ++col)
00642     {
00643           eT* t_coldata = t.slice_colptr(t.aux_slice1, col);
00644     const eT* x_coldata = x.colptr(col);
00645     
00646     for(u32 row = 0; row < t.n_rows; ++row)
00647       {
00648       t_coldata[row] -= x_coldata[row];
00649       }
00650     }
00651   }
00652 
00653 
00654 
00655 template<typename eT>
00656 template<typename T1>
00657 inline
00658 void
00659 subview_cube<eT>::operator%= (const Base<eT,T1>& in)
00660   {
00661   arma_extra_debug_sigprint();
00662   
00663   const unwrap<T1> tmp(in.get_ref());
00664   
00665   const Mat<eT>&          x = tmp.M;
00666         subview_cube<eT>& t = *this;
00667   
00668   arma_debug_assert_same_size(t, x, "cube schur product");
00669   
00670   for(u32 col = 0; col<t.n_cols; ++col)
00671     {
00672           eT* t_coldata = t.slice_colptr(t.aux_slice1, col);
00673     const eT* x_coldata = x.colptr(col);
00674     
00675     for(u32 row = 0; row<t.n_rows; ++row)
00676       {
00677       t_coldata[row] *= x_coldata[row];
00678       }
00679     }
00680   }
00681 
00682 
00683 
00684 template<typename eT>
00685 template<typename T1>
00686 inline
00687 void
00688 subview_cube<eT>::operator/= (const Base<eT,T1>& in)
00689   {
00690   arma_extra_debug_sigprint();
00691   
00692   const unwrap<T1> tmp(in.get_ref());
00693   
00694   const Mat<eT>&          x = tmp.M;
00695         subview_cube<eT>& t = *this;
00696   
00697   arma_debug_assert_same_size(t, x, "element-wise cube division");
00698   
00699   
00700   for(u32 col = 0; col<t.n_cols; ++col)
00701     {
00702           eT* t_coldata = t.slice_colptr(t.aux_slice1, col);
00703     const eT* x_coldata = x.colptr(col);
00704     
00705     for(u32 row = 0; row<t.n_rows; ++row)
00706       {
00707       t_coldata[row] /= x_coldata[row];
00708       }
00709     }
00710   }
00711 
00712 
00713 
00714 template<typename eT>
00715 inline
00716 void
00717 subview_cube<eT>::fill(const eT val)
00718   {
00719   arma_extra_debug_sigprint();
00720 
00721   for(u32 slice = 0; slice < n_slices; ++slice)
00722     {
00723     for(u32 col = 0; col < n_cols; ++col)
00724       {
00725       eT* coldata = slice_colptr(slice,col);
00726       
00727       for(u32 row = 0; row < n_rows; ++row)
00728         {
00729         coldata[row] = val;
00730         }
00731       }
00732     }
00733   
00734   }
00735 
00736 
00737 
00738 template<typename eT>
00739 inline
00740 void
00741 subview_cube<eT>::zeros()
00742   {
00743   arma_extra_debug_sigprint();
00744   
00745   fill(eT(0));
00746   }
00747 
00748 
00749 
00750 template<typename eT>
00751 inline
00752 void
00753 subview_cube<eT>::ones()
00754   {
00755   arma_extra_debug_sigprint();
00756   
00757   fill(eT(1));
00758   }
00759 
00760 
00761 
00762 template<typename eT>
00763 arma_inline
00764 eT&
00765 subview_cube<eT>::operator[](const u32 i)
00766   {
00767   arma_check( (m_ptr == 0), "subview_cube::operator[]: cube is read-only");
00768   
00769   const u32 in_slice = i / n_elem_slice;
00770   const u32 offset   = in_slice * n_elem_slice;
00771   const u32 j        = i - offset;
00772   
00773   const u32 in_col   = j / n_rows;
00774   const u32 in_row   = j % n_rows;
00775 
00776   const u32 index = (in_slice + aux_slice1)*m.n_elem_slice + (in_col + aux_col1)*m.n_rows + aux_row1 + in_row;
00777   return access::rw( (*m_ptr).mem[index] );
00778   }
00779 
00780 
00781 
00782 template<typename eT>
00783 arma_inline
00784 eT
00785 subview_cube<eT>::operator[](const u32 i) const
00786   {
00787   const u32 in_slice = i / n_elem_slice;
00788   const u32 offset   = in_slice * n_elem_slice;
00789   const u32 j        = i - offset;
00790   
00791   const u32 in_col   = j / n_rows;
00792   const u32 in_row   = j % n_rows;
00793 
00794   const u32 index = (in_slice + aux_slice1)*m.n_elem_slice + (in_col + aux_col1)*m.n_rows + aux_row1 + in_row;
00795   return m.mem[index];
00796   }
00797 
00798 
00799 
00800 template<typename eT>
00801 arma_inline
00802 eT&
00803 subview_cube<eT>::operator()(const u32 i)
00804   {
00805   arma_check( (m_ptr == 0), "subview_cube::operator(): matrix is read-only");
00806   arma_debug_check( (i >= n_elem), "subview_cube::operator(): index out of bounds");
00807   
00808   const u32 in_slice = i / n_elem_slice;
00809   const u32 offset   = in_slice * n_elem_slice;
00810   const u32 j        = i - offset;
00811   
00812   const u32 in_col   = j / n_rows;
00813   const u32 in_row   = j % n_rows;
00814 
00815   const u32 index = (in_slice + aux_slice1)*m.n_elem_slice + (in_col + aux_col1)*m.n_rows + aux_row1 + in_row;
00816   return access::rw( (*m_ptr).mem[index] );
00817   }
00818 
00819 
00820 
00821 template<typename eT>
00822 arma_inline
00823 eT
00824 subview_cube<eT>::operator()(const u32 i) const
00825   {
00826   arma_debug_check( (i >= n_elem), "subview_cube::operator(): index out of bounds");
00827   
00828   const u32 in_slice = i / n_elem_slice;
00829   const u32 offset   = in_slice * n_elem_slice;
00830   const u32 j        = i - offset;
00831   
00832   const u32 in_col   = j / n_rows;
00833   const u32 in_row   = j % n_rows;
00834 
00835   const u32 index = (in_slice + aux_slice1)*m.n_elem_slice + (in_col + aux_col1)*m.n_rows + aux_row1 + in_row;
00836   return m.mem[index];
00837   }
00838 
00839 
00840 
00841 template<typename eT>
00842 arma_inline
00843 eT&
00844 subview_cube<eT>::operator()(const u32 in_row, const u32 in_col, const u32 in_slice)
00845   {
00846   arma_check( (m_ptr == 0), "subview_cube::operator(): matrix is read-only");
00847   arma_debug_check( ( (in_row >= n_rows) || (in_col >= n_cols) || (in_slice >= n_slices) ), "subview_cube::operator(): location out of bounds");
00848   
00849   const u32 index = (in_slice + aux_slice1)*m.n_elem_slice + (in_col + aux_col1)*m.n_rows + aux_row1 + in_row;
00850   return access::rw( (*m_ptr).mem[index] );
00851   }
00852 
00853 
00854 
00855 template<typename eT>
00856 arma_inline
00857 eT
00858 subview_cube<eT>::operator()(const u32 in_row, const u32 in_col, const u32 in_slice) const
00859   {
00860   arma_debug_check( ( (in_row >= n_rows) || (in_col >= n_cols) || (in_slice >= n_slices) ), "subview_cube::operator(): location out of bounds");
00861   
00862   const u32 index = (in_slice + aux_slice1)*m.n_elem_slice + (in_col + aux_col1)*m.n_rows + aux_row1 + in_row;
00863   return m.mem[index];
00864   }
00865 
00866 
00867 
00868 template<typename eT>
00869 arma_inline
00870 eT&
00871 subview_cube<eT>::at(const u32 in_row, const u32 in_col, const u32 in_slice)
00872   {
00873   arma_check( (m_ptr == 0), "subview_cube::at(): cube is read-only");
00874   
00875   const u32 index = (in_slice + aux_slice1)*m.n_elem_slice + (in_col + aux_col1)*m.n_rows + aux_row1 + in_row;
00876   return access::rw( (*m_ptr).mem[index] );
00877   }
00878 
00879 
00880 
00881 template<typename eT>
00882 arma_inline
00883 eT
00884 subview_cube<eT>::at(const u32 in_row, const u32 in_col, const u32 in_slice) const
00885   {
00886   const u32 index = (in_slice + aux_slice1)*m.n_elem_slice + (in_col + aux_col1)*m.n_rows + aux_row1 + in_row;
00887   return m.mem[index];
00888   }
00889 
00890 
00891 
00892 template<typename eT>
00893 arma_inline
00894 eT*
00895 subview_cube<eT>::slice_colptr(const u32 in_slice, const u32 in_col)
00896   {
00897   arma_check( (m_ptr == 0), "subview_cube::slice_colptr(): cube is read-only");
00898     
00899   return & access::rw((*m_ptr).mem[  (in_slice + aux_slice1)*m.n_elem_slice + (in_col + aux_col1)*m.n_rows + aux_row1  ]);
00900   }
00901 
00902 
00903 
00904 template<typename eT>
00905 arma_inline
00906 const eT*
00907 subview_cube<eT>::slice_colptr(const u32 in_slice, const u32 in_col) const
00908   {
00909   return & m.mem[ (in_slice + aux_slice1)*m.n_elem_slice + (in_col + aux_col1)*m.n_rows + aux_row1 ];
00910   }
00911 
00912 
00913 
00914 template<typename eT>
00915 inline
00916 bool
00917 subview_cube<eT>::check_overlap(const subview_cube<eT>& x) const
00918   {
00919   const subview_cube<eT>& t = *this;
00920   
00921   if(&t.m != &x.m)
00922     {
00923     return false;
00924     }
00925   else
00926     {
00927     const bool row_overlap =
00928       (
00929       ( (x.aux_row1 >= t.aux_row1) && (x.aux_row1 <= t.aux_row2) )
00930       || 
00931       ( (x.aux_row2 >= t.aux_row1) && (x.aux_row2 <= t.aux_row2) )
00932       );
00933     
00934     const bool col_overlap =
00935       (
00936       ( (x.aux_col1 >= t.aux_col1) && (x.aux_col1 <= t.aux_col2) )
00937       || 
00938       ( (x.aux_col2 >= t.aux_col1) && (x.aux_col2 <= t.aux_col2) )
00939       );
00940     
00941     const bool slice_overlap =
00942       (
00943       ( (x.aux_slice1 >= t.aux_slice1) && (x.aux_slice1 <= t.aux_slice2) )
00944       || 
00945       ( (x.aux_slice2 >= t.aux_slice1) && (x.aux_slice2 <= t.aux_slice2) )
00946       );
00947     
00948     return (row_overlap & col_overlap & slice_overlap);
00949     }
00950   }
00951 
00952 
00953 
00954 template<typename eT>
00955 inline
00956 bool
00957 subview_cube<eT>::check_overlap(const Mat<eT>& x) const
00958   {
00959   const subview_cube<eT>& t = *this;
00960   
00961   for(u32 slice = t.aux_slice1; slice <= t.aux_slice2; ++slice)
00962     {
00963     const Mat<eT>& y = *(t.m.mat_ptrs[slice]);
00964   
00965     if( x.memptr() == y.memptr() )
00966       {
00967       return true;
00968       }
00969     }
00970   
00971   return false;
00972   }
00973 
00974 
00975 
00976 //! cube X = Y.subcube(...)
00977 template<typename eT>
00978 inline
00979 void
00980 subview_cube<eT>::extract(Cube<eT>& actual_out, const subview_cube<eT>& in)
00981   {
00982   arma_extra_debug_sigprint();
00983   
00984   //
00985   const bool alias = (&actual_out == &in.m);
00986   
00987   Cube<eT>* tmp = (alias) ? new Cube<eT> : 0;
00988   Cube<eT>& out = (alias) ? (*tmp)       : actual_out;
00989   
00990   //
00991   
00992   const u32 n_rows   = in.n_rows;
00993   const u32 n_cols   = in.n_cols;
00994   const u32 n_slices = in.n_slices;
00995   
00996   out.set_size(n_rows, n_cols, n_slices);
00997   
00998   arma_extra_debug_print(arma_boost::format("out.n_rows = %d   out.n_cols = %d    out.n_slices = %d    in.m.n_rows = %d   in.m.n_cols = %d   in.m.n_slices = %d") % out.n_rows % out.n_cols % out.n_slices % in.m.n_rows % in.m.n_cols % in.m.n_slices);
00999   
01000   
01001   for(u32 slice = 0; slice<n_slices; ++slice)
01002     {
01003     for(u32 col = 0; col<n_cols; ++col)
01004       {
01005             eT* out_coldata = out.slice_colptr(slice,col);
01006       const eT*  in_coldata =  in.slice_colptr(slice,col);
01007       
01008       for(u32 row = 0; row<n_rows; ++row)
01009         {
01010         out_coldata[row] = in_coldata[row];
01011         }
01012       
01013       }
01014     }
01015   
01016   
01017   if(alias)
01018     {
01019     actual_out = out;
01020     delete tmp;
01021     }
01022   
01023   }
01024 
01025 
01026 
01027 //! mat X = Y.subcube(...)
01028 template<typename eT>
01029 inline
01030 void
01031 subview_cube<eT>::extract(Mat<eT>& out, const subview_cube<eT>& in)
01032   {
01033   arma_extra_debug_sigprint();
01034   
01035   arma_debug_check( (in.n_slices != 1), "subview_cube::extract(): given subcube doesn't have exactly one slice" );
01036   
01037   out.set_size(in.n_rows, in.n_cols);
01038   
01039   for(u32 col = 0; col < in.n_cols; ++col)
01040     {
01041     const eT* in_coldata  = in.slice_colptr(in.aux_slice1, col);
01042           eT* out_coldata = out.colptr(col);
01043     
01044     for(u32 row = 0; row < in.n_rows; ++row)
01045       {
01046       out_coldata[row] = in_coldata[row];
01047       }
01048     }
01049   }
01050 
01051 
01052 
01053 //! cube X += Y.subcube(...)
01054 template<typename eT>
01055 inline
01056 void
01057 subview_cube<eT>::plus_inplace(Cube<eT>& out, const subview_cube<eT>& in)
01058   {
01059   arma_extra_debug_sigprint();
01060   
01061   arma_debug_assert_same_size(out, in, "cube addition");
01062   
01063   if(&out != &in.m)
01064     {
01065     const u32 n_rows   = out.n_rows;
01066     const u32 n_cols   = out.n_cols;
01067     const u32 n_slices = out.n_slices;
01068     
01069     for(u32 slice = 0; slice<n_slices; ++slice)
01070       {
01071       for(u32 col = 0; col<n_cols; ++col)
01072         {
01073               eT* out_coldata = out.slice_colptr(slice,col);
01074         const eT*  in_coldata =  in.slice_colptr(slice,col);
01075         
01076         for(u32 row = 0; row<n_rows; ++row)
01077           {
01078           out_coldata[row] += in_coldata[row];
01079           }
01080         }
01081       }
01082     }
01083   else
01084     {
01085     // X += X.subcube(...)
01086     // this only makes sense if X and X.subcube(...) are the same size
01087     
01088     eT* out_mem = out.memptr();
01089     
01090     for(u32 i=0; i<out.n_elem; ++i)
01091       {
01092       out_mem[i] *= eT(2);
01093       }
01094     }
01095   
01096   }
01097 
01098 
01099 
01100 //! cube X -= Y.subcube(...)
01101 template<typename eT>
01102 inline
01103 void
01104 subview_cube<eT>::minus_inplace(Cube<eT>& out, const subview_cube<eT>& in)
01105   {
01106   arma_extra_debug_sigprint();
01107   
01108   arma_debug_assert_same_size(out, in, "cube subtraction");
01109   
01110   if(&out != &in.m)
01111     {
01112     const u32 n_rows   = out.n_rows;
01113     const u32 n_cols   = out.n_cols;
01114     const u32 n_slices = out.n_slices;
01115     
01116     for(u32 slice = 0; slice<n_slices; ++slice)
01117       {
01118       for(u32 col = 0; col<n_cols; ++col)
01119         {
01120               eT* out_coldata = out.slice_colptr(slice,col);
01121         const eT*  in_coldata =  in.slice_colptr(slice,col);
01122         
01123         for(u32 row = 0; row<n_rows; ++row)
01124           {
01125           out_coldata[row] -= in_coldata[row];
01126           }
01127         }
01128       }
01129     }
01130   else
01131     {
01132     // X -= X.subcube(...)
01133     // this only makes sense if X and X.subcube(...) are the same size
01134 
01135     out.zeros();
01136     }
01137   
01138   }
01139 
01140 
01141 
01142 //! cube X %= Y.subcube(...)
01143 template<typename eT>
01144 inline
01145 void
01146 subview_cube<eT>::schur_inplace(Cube<eT>& out, const subview_cube<eT>& in)
01147   {
01148   arma_extra_debug_sigprint();
01149   
01150   arma_debug_assert_same_size(out, in, "cube schur product");
01151   
01152   if(&out != &in.m)
01153     {
01154     const u32 n_rows   = out.n_rows;
01155     const u32 n_cols   = out.n_cols;
01156     const u32 n_slices = out.n_slices;
01157     
01158     for(u32 slice = 0; slice<n_slices; ++slice)
01159       {
01160       for(u32 col = 0; col<n_cols; ++col)
01161         {
01162               eT* out_coldata = out.slice_colptr(slice,col);
01163         const eT*  in_coldata =  in.slice_colptr(slice,col);
01164         
01165         for(u32 row = 0; row<n_rows; ++row)
01166           {
01167           out_coldata[row] *= in_coldata[row];
01168           }
01169         }
01170       }
01171     }
01172   else
01173     {
01174     // X %= X.subcube(...)
01175     // this only makes sense if X and X.subcube(...) are the same size
01176 
01177     eT* out_mem = out.memptr();
01178     
01179     for(u32 i=0; i<out.n_elem; ++i)
01180       {
01181       const eT tmp = out_mem[i];
01182       out_mem[i] = tmp*tmp;
01183       }
01184     }
01185   
01186   }
01187 
01188 
01189 
01190 //! cube X /= Y.subcube(...)
01191 template<typename eT>
01192 inline
01193 void
01194 subview_cube<eT>::div_inplace(Cube<eT>& out, const subview_cube<eT>& in)
01195   {
01196   arma_extra_debug_sigprint();
01197   
01198   arma_debug_assert_same_size(out, in, "element-wise cube division");
01199   
01200   if(&out != &in.m)
01201     {
01202     const u32 n_rows   = out.n_rows;
01203     const u32 n_cols   = out.n_cols;
01204     const u32 n_slices = out.n_slices;
01205     
01206     for(u32 slice = 0; slice<n_slices; ++slice)
01207       {
01208       for(u32 col = 0; col<n_cols; ++col)
01209         {
01210               eT* out_coldata = out.slice_colptr(slice,col);
01211         const eT*  in_coldata =  in.slice_colptr(slice,col);
01212         
01213         for(u32 row = 0; row<n_rows; ++row)
01214           {
01215           out_coldata[row] /= in_coldata[row];
01216           }
01217         }
01218       }
01219     }
01220   else
01221     {
01222     // X /= X.subcube(...)
01223     // this only makes sense if X and X.subcube(...) are the same size
01224 
01225     eT* out_mem = out.memptr();
01226     
01227     for(u32 i=0; i<out.n_elem; ++i)
01228       {
01229       const eT tmp = out_mem[i];
01230       out_mem[i] = tmp/tmp;  // using tmp/tmp as tmp might be zero
01231       }
01232     }
01233   
01234   }
01235 
01236 
01237 
01238 //! mat X += Y.subcube(...)
01239 template<typename eT>
01240 inline
01241 void
01242 subview_cube<eT>::plus_inplace(Mat<eT>& out, const subview_cube<eT>& in)
01243   {
01244   arma_extra_debug_sigprint();
01245   
01246   arma_debug_assert_same_size(out, in, "matrix addition");
01247   
01248   for(u32 col = 0; col < in.n_cols; ++col)
01249     {
01250     const eT* in_coldata  = in.slice_colptr(in.aux_slice1, col);
01251           eT* out_coldata = out.colptr(col);
01252     
01253     for(u32 row = 0; row < in.n_rows; ++row)
01254       {
01255       out_coldata[row] += in_coldata[row];
01256       }
01257     }
01258     
01259   }
01260 
01261 
01262 
01263 //! mat X -= Y.subcube(...)
01264 template<typename eT>
01265 inline
01266 void
01267 subview_cube<eT>::minus_inplace(Mat<eT>& out, const subview_cube<eT>& in)
01268   {
01269   arma_extra_debug_sigprint();
01270   
01271   arma_debug_assert_same_size(out, in, "matrix subtraction");
01272   
01273   for(u32 col = 0; col < in.n_cols; ++col)
01274     {
01275     const eT* in_coldata  = in.slice_colptr(in.aux_slice1, col);
01276           eT* out_coldata = out.colptr(col);
01277     
01278     for(u32 row = 0; row < in.n_rows; ++row)
01279       {
01280       out_coldata[row] -= in_coldata[row];
01281       }
01282     }
01283     
01284   }
01285 
01286 
01287 
01288 //! mat X %= Y.subcube(...)
01289 template<typename eT>
01290 inline
01291 void
01292 subview_cube<eT>::schur_inplace(Mat<eT>& out, const subview_cube<eT>& in)
01293   {
01294   arma_extra_debug_sigprint();
01295   
01296   arma_debug_assert_same_size(out, in, "matrix schur product");
01297   
01298   for(u32 col = 0; col < in.n_cols; ++col)
01299     {
01300     const eT* in_coldata  = in.slice_colptr(in.aux_slice1, col);
01301           eT* out_coldata = out.colptr(col);
01302     
01303     for(u32 row = 0; row < in.n_rows; ++row)
01304       {
01305       out_coldata[row] *= in_coldata[row];
01306       }
01307     }
01308     
01309   }
01310 
01311 
01312 
01313 //! mat X /= Y.subcube(...)
01314 template<typename eT>
01315 inline
01316 void
01317 subview_cube<eT>::div_inplace(Mat<eT>& out, const subview_cube<eT>& in)
01318   {
01319   arma_extra_debug_sigprint();
01320   
01321   arma_debug_assert_same_size(out, in, "matrix element-wise division");
01322   
01323   for(u32 col = 0; col < in.n_cols; ++col)
01324     {
01325     const eT* in_coldata  = in.slice_colptr(in.aux_slice1, col);
01326           eT* out_coldata = out.colptr(col);
01327     
01328     for(u32 row = 0; row < in.n_rows; ++row)
01329       {
01330       out_coldata[row] /= in_coldata[row];
01331       }
01332     }
01333     
01334   }
01335 
01336 
01337 
01338 //! @}