1 /*
2  *  position.h
3  *
4  *  This file is part of NEST.
5  *
6  *  Copyright (C) 2004 The NEST Initiative
7  *
8  *  NEST is free software: you can redistribute it and/or modify
9  *  it under the terms of the GNU General Public License as published by
10  *  the Free Software Foundation, either version 2 of the License, or
11  *  (at your option) any later version.
12  *
13  *  NEST is distributed in the hope that it will be useful,
14  *  but WITHOUT ANY WARRANTY; without even the implied warranty of
15  *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
16  *  GNU General Public License for more details.
17  *
18  *  You should have received a copy of the GNU General Public License
19  *  along with NEST.  If not, see <http://www.gnu.org/licenses/>.
20  *
21  */
22 
23 #ifndef POSITION_H
24 #define POSITION_H
25 
26 // C++ includes:
27 #include <algorithm>
28 #include <array>
29 #include <cassert>
30 #include <cmath>
31 #include <iostream>
32 #include <sstream>
33 #include <string>
34 #include <utility>
35 #include <vector>
36 
37 // Includes from libnestutil:
38 #include "compose.hpp"
39 
40 // Includes from nestkernel:
41 #include "exceptions.h"
42 #include "nest_types.h"
43 
44 // Includes from sli:
45 #include "token.h"
46 
47 namespace nest
48 {
49 
50 // It is necessary to declare the template for operator<< first in order
51 // to get the friend declaration to work
52 template < int D, class T >
53 class Position;
54 
55 template < int D, class T >
56 std::ostream& operator<<( std::ostream& os, const Position< D, T >& pos );
57 
58 template < int D, class T = double >
59 class Position
60 {
61 public:
62   template < int OD, class OT >
63   friend class Position;
64 
65   /**
66    * Default constructor, initializing all coordinates to zero.
67    */
68   Position();
69 
70   /**
71    * 2D Constructor.
72    */
73   Position( const T&, const T& );
74 
75   /**
76    * 3D Constructor.
77    */
78   Position( const T&, const T&, const T& );
79 
80   /**
81    * Constructor initializing a Position from an array.
82    */
83   Position( const T* const y );
84 
85   /**
86    * Constructor initializing a Position from a std::vector.
87    */
88   Position( const std::vector< T >& y );
89 
90   /**
91    * Copy constructor.
92    */
93   Position( const Position& other );
94 
95   template < class U >
96   Position( const Position< D, U >& other );
97 
98   /**
99    * Move constructor.
100    */
101   Position( Position&& other );
102 
103   /**
104    * Assignment constructor.
105    */
106   Position& operator=( const Position& other ) = default;
107   Position& operator=( const std::vector< T >& y );
108 
109   /**
110    * Move assignment constructor.
111    */
112   Position& operator=( Position&& other ) = default;
113 
114   /**
115    * @returns an element (coordinate) of the Position
116    */
117   T& operator[]( int i );
118 
119   /**
120    * @returns an element (coordinate) of the Position
121    */
122   const T& operator[]( int i ) const;
123 
124   /**
125    * Moves Position variables into an array.
126    * @returns array of positions stored as a token object.
127    */
128   Token getToken() const;
129 
130   const std::vector< T > get_vector() const;
131   void get_vector( std::vector< T >& vector ) const;
132 
133   /**
134    * Elementwise addition.
135    * @returns elementwise sum of coordinates.
136    */
137   template < class OT >
138   Position operator+( const Position< D, OT >& other ) const;
139 
140   /**
141    * Elementwise subtraction.
142    * @returns elementwise difference of coordinates.
143    */
144   template < class OT >
145   Position operator-( const Position< D, OT >& other ) const;
146 
147   /**
148    * Unary minus.
149    * @returns opposite vector.
150    */
151   Position operator-() const;
152 
153   /**
154    * Elementwise multiplication.
155    * @returns elementwise product of coordinates.
156    */
157   template < class OT >
158   Position operator*( const Position< D, OT >& other ) const;
159 
160   /**
161    * Elementwise division.
162    * @returns elementwise quotient of coordinates.
163    */
164   template < class OT >
165   Position operator/( const Position< D, OT >& other ) const;
166 
167   /**
168    * Elementwise addition with scalar
169    * @returns position vector with scalar added to all coordinates
170    */
171   Position operator+( const T& ) const;
172 
173   /**
174    * Elementwise subtraction with scalar
175    * @returns position vector with scalar subtracted from all coordinates
176    */
177   Position operator-( const T& ) const;
178 
179   /**
180    * Multiplication with scalar
181    * @returns position vector multiplied with the scalar.
182    */
183   Position operator*( const T& ) const;
184 
185   /**
186    * Division with scalar
187    * @returns position vector divided by the scalar.
188    */
189   Position operator/( const T& ) const;
190 
191   /**
192    * In-place elementwise addition.
193    * @returns the Position itself after adding the other Position
194    * elementwise.
195    */
196   template < class OT >
197   Position& operator+=( const Position< D, OT >& );
198 
199   /**
200    * In-place elementwise subtraction.
201    * @returns the Position itself after subtracting the other Position
202    * elementwise.
203    */
204   template < class OT >
205   Position& operator-=( const Position< D, OT >& );
206 
207   /**
208    * In-place elementwise multiplication.
209    * @returns the Position itself after multiplying with the other
210    * Position elementwise.
211    */
212   template < class OT >
213   Position& operator*=( const Position< D, OT >& );
214 
215   /**
216    * In-place elementwise division.
217    * @returns the Position itself after dividing by the other Position
218    * elementwise.
219    */
220   template < class OT >
221   Position& operator/=( const Position< D, OT >& );
222 
223   /**
224    * In-place elementwise addition with scalar.
225    * @returns the Position itself after adding the scalar to all coordinates.
226    */
227   Position& operator+=( const T& );
228 
229   /**
230    * In-place elementwise subtraction with scalar.
231    * @returns the Position itself after subtracting the scalar from all
232    * coordinates.
233    */
234   Position& operator-=( const T& );
235 
236   /**
237    * In-place multiplication by scalar.
238    * @returns the Position itself after multiplying with the scalar.
239    */
240   Position& operator*=( const T& );
241 
242   /**
243    * In-place elementwise division.
244    * @returns the Position itself after dividing by the scalar.
245    */
246   Position& operator/=( const T& );
247 
248   /**
249    * @returns true if all coordinates are equal
250    */
251   bool operator==( const Position& y ) const;
252 
253   /**
254    * @returns true if not all coordinates are equal
255    */
256   bool operator!=( const Position& y ) const;
257 
258   /**
259    * @returns true if all coordinates are less
260    */
261   bool operator<( const Position& y ) const;
262 
263   /**
264    * @returns true if all coordinates are greater
265    */
266   bool operator>( const Position& y ) const;
267 
268   /**
269    * @returns true if all coordinates are less or equal
270    */
271   bool operator<=( const Position& y ) const;
272 
273   /**
274    * @returns true if all coordinates are greater or equal
275    */
276   bool operator>=( const Position& y ) const;
277 
278   /**
279    * Length of Position vector.
280    * @returns Euclidian norm of the vector.
281    */
282   T length() const;
283 
284   /**
285    * @returns string representation of Position
286    */
287   operator std::string() const;
288 
289   /**
290    * Print position to output stream.
291    *
292    * Format: Only as many coordinates as dimensions,
293    *         separated by spaces [default], no trailing space.
294    *
295    * @param out output stream
296    * @param sep separator character
297    */
298   void print( std::ostream& out, char sep = ' ' ) const;
299 
300   /**
301    * Output the Position to an ostream.
302    */
303   friend std::ostream& operator<<<>( std::ostream& os, const Position< D, T >& pos );
304 
305 protected:
306   std::array< T, D > x_;
307 };
308 
309 /**
310  * A box is defined by the lower left corner (minimum coordinates) and
311  * the upper right corner (maximum coordinates).
312  */
313 template < int D >
314 struct Box
315 {
BoxBox316   Box()
317   {
318   }
BoxBox319   Box( const Position< D >& lower_left, const Position< D >& upper_right )
320     : lower_left( lower_left )
321     , upper_right( upper_right )
322   {
323   }
324 
325   Position< D > lower_left;
326   Position< D > upper_right;
327 };
328 
329 /**
330  * An index into a multidimensional array.
331  */
332 template < int D >
333 class MultiIndex : public Position< D, int >
334 {
335 public:
MultiIndex()336   MultiIndex()
337     : Position< D, int >()
338     , lower_left_()
339     , upper_right_()
340   {
341   }
342 
MultiIndex(const Position<D,int> & ur)343   MultiIndex( const Position< D, int >& ur )
344     : Position< D, int >()
345     , lower_left_()
346     , upper_right_( ur )
347   {
348   }
349 
MultiIndex(const Position<D,int> & lower_left,const Position<D,int> & upper_right)350   MultiIndex( const Position< D, int >& lower_left, const Position< D, int >& upper_right )
351     : Position< D, int >( lower_left )
352     , lower_left_( lower_left )
353     , upper_right_( upper_right )
354   {
355   }
356 
357   MultiIndex& operator++()
358   {
359     // Try increasing the first coordinate first, resetting it and
360     // continuing with the next if the first one overflows, and so on
361     for ( int i = 0; i < D; ++i )
362     {
363       this->x_[ i ]++;
364       if ( this->x_[ i ] < upper_right_[ i ] )
365       {
366         return *this;
367       }
368       this->x_[ i ] = lower_left_[ i ];
369     }
370     // If we reach this point, we are outside of bounds. The upper
371     // right point is used as a marker to show that we have reached the
372     // end.
373     for ( int i = 0; i < D; ++i )
374     {
375       this->x_[ i ] = upper_right_[ i ];
376     }
377     return *this;
378   }
379 
380   MultiIndex operator++( int )
381   {
382     MultiIndex tmp = *this;
383     ++*this;
384     return tmp;
385   }
386 
387   Position< D, int >
get_lower_left()388   get_lower_left() const
389   {
390     return lower_left_;
391   }
392   Position< D, int >
get_upper_right()393   get_upper_right() const
394   {
395     return upper_right_;
396   }
397 
398 private:
399   Position< D, int > lower_left_;
400   Position< D, int > upper_right_;
401 };
402 
403 
404 template < int D, class T >
Position()405 inline Position< D, T >::Position()
406 {
407   x_.fill( 0 );
408 }
409 
410 template < int D, class T >
Position(const T & x,const T & y)411 inline Position< D, T >::Position( const T& x, const T& y )
412 {
413   assert( D == 2 );
414   x_[ 0 ] = x;
415   x_[ 1 ] = y;
416 }
417 
418 template < int D, class T >
Position(const T & x,const T & y,const T & z)419 inline Position< D, T >::Position( const T& x, const T& y, const T& z )
420 {
421   assert( D == 3 );
422   x_[ 0 ] = x;
423   x_[ 1 ] = y;
424   x_[ 2 ] = z;
425 }
426 
427 template < int D, class T >
Position(const T * const y)428 inline Position< D, T >::Position( const T* const y )
429 {
430   for ( int i = 0; i < D; ++i )
431   {
432     x_[ i ] = y[ i ];
433   }
434 }
435 
436 template < int D, class T >
Position(const std::vector<T> & y)437 inline Position< D, T >::Position( const std::vector< T >& y )
438 {
439   if ( y.size() != D )
440   {
441     throw BadProperty( String::compose( "Expected a %1-dimensional position.", D ) );
442   }
443   std::copy( y.begin(), y.end(), x_.begin() );
444 }
445 
446 template < int D, class T >
Position(const Position<D,T> & other)447 inline Position< D, T >::Position( const Position< D, T >& other )
448   : x_( other.x_ )
449 {
450 }
451 
452 template < int D, class T >
453 template < class U >
Position(const Position<D,U> & other)454 inline Position< D, T >::Position( const Position< D, U >& other )
455   : x_( other.x_ )
456 {
457 }
458 
459 template < int D, class T >
Position(Position && other)460 inline Position< D, T >::Position( Position&& other )
461 {
462   x_ = std::move( other.x_ );
463 }
464 
465 template < int D, class T >
466 inline Position< D, T >& Position< D, T >::operator=( const std::vector< T >& y )
467 {
468   if ( y.size() != D )
469   {
470     throw BadProperty( String::compose( "Expected a %1-dimensional position.", D ) );
471   }
472   std::copy( y.begin(), y.end(), x_.begin() );
473 
474   return *this;
475 }
476 
477 template < int D, class T >
478 inline T& Position< D, T >::operator[]( int i )
479 {
480   return x_[ i ];
481 }
482 
483 template < int D, class T >
484 inline const T& Position< D, T >::operator[]( int i ) const
485 {
486   return x_[ i ];
487 }
488 
489 template < int D, class T >
490 Token
getToken()491 Position< D, T >::getToken() const
492 {
493   std::vector< T > result = get_vector();
494   return Token( result );
495 }
496 
497 
498 template < int D, class T >
499 const std::vector< T >
get_vector()500 Position< D, T >::get_vector() const
501 {
502   return std::vector< T >( x_.begin(), x_.end() );
503 }
504 
505 template < int D, class T >
506 void
get_vector(std::vector<T> & vector)507 Position< D, T >::get_vector( std::vector< T >& vector ) const
508 {
509   assert( vector.size() == D );
510   std::copy( x_.begin(), x_.end(), vector.begin() );
511 }
512 
513 
514 template < int D, class T >
515 template < class OT >
516 inline Position< D, T > Position< D, T >::operator+( const Position< D, OT >& other ) const
517 {
518   Position p = *this;
519   p += other;
520   return p;
521 }
522 
523 template < int D, class T >
524 template < class OT >
525 inline Position< D, T > Position< D, T >::operator-( const Position< D, OT >& other ) const
526 {
527   Position p = *this;
528   p -= other;
529   return p;
530 }
531 
532 template < int D, class T >
533 inline Position< D, T > Position< D, T >::operator-() const
534 {
535   Position p;
536   p -= *this;
537   return p;
538 }
539 
540 template < int D, class T >
541 template < class OT >
542 inline Position< D, T > Position< D, T >::operator*( const Position< D, OT >& other ) const
543 {
544   Position p = *this;
545   p *= other;
546   return p;
547 }
548 
549 template < int D, class T >
550 template < class OT >
551 inline Position< D, T > Position< D, T >::operator/( const Position< D, OT >& other ) const
552 {
553   Position p = *this;
554   p /= other;
555   return p;
556 }
557 
558 template < int D, class T >
559 inline Position< D, T > Position< D, T >::operator+( const T& a ) const
560 {
561   Position p = *this;
562   p += a;
563   return p;
564 }
565 
566 template < int D, class T >
567 inline Position< D, T > Position< D, T >::operator-( const T& a ) const
568 {
569   Position p = *this;
570   p -= a;
571   return p;
572 }
573 
574 template < int D, class T >
575 inline Position< D, T > Position< D, T >::operator*( const T& a ) const
576 {
577   Position p = *this;
578   p *= a;
579   return p;
580 }
581 
582 template < int D, class T >
583 inline Position< D, T > Position< D, T >::operator/( const T& a ) const
584 {
585   Position p = *this;
586   p /= a;
587   return p;
588 }
589 
590 template < int D, class T >
591 template < class OT >
592 inline Position< D, T >& Position< D, T >::operator+=( const Position< D, OT >& other )
593 {
594   for ( int i = 0; i < D; ++i )
595   {
596     x_[ i ] += other.x_[ i ];
597   }
598   return *this;
599 }
600 
601 template < int D, class T >
602 template < class OT >
603 inline Position< D, T >& Position< D, T >::operator-=( const Position< D, OT >& other )
604 {
605   for ( int i = 0; i < D; ++i )
606   {
607     x_[ i ] -= other.x_[ i ];
608   }
609   return *this;
610 }
611 
612 template < int D, class T >
613 template < class OT >
614 inline Position< D, T >& Position< D, T >::operator*=( const Position< D, OT >& other )
615 {
616   for ( int i = 0; i < D; ++i )
617   {
618     x_[ i ] *= other.x_[ i ];
619   }
620   return *this;
621 }
622 
623 template < int D, class T >
624 template < class OT >
625 inline Position< D, T >& Position< D, T >::operator/=( const Position< D, OT >& other )
626 {
627   for ( int i = 0; i < D; ++i )
628   {
629     x_[ i ] /= other.x_[ i ];
630   }
631   return *this;
632 }
633 
634 template < int D, class T >
635 inline Position< D, T >& Position< D, T >::operator+=( const T& a )
636 {
637   for ( int i = 0; i < D; ++i )
638   {
639     x_[ i ] += a;
640   }
641   return *this;
642 }
643 
644 template < int D, class T >
645 inline Position< D, T >& Position< D, T >::operator-=( const T& a )
646 {
647   for ( int i = 0; i < D; ++i )
648   {
649     x_[ i ] -= a;
650   }
651   return *this;
652 }
653 
654 template < int D, class T >
655 inline Position< D, T >& Position< D, T >::operator*=( const T& a )
656 {
657   for ( int i = 0; i < D; ++i )
658   {
659     x_[ i ] *= a;
660   }
661   return *this;
662 }
663 
664 template < int D, class T >
665 inline Position< D, T >& Position< D, T >::operator/=( const T& a )
666 {
667   for ( int i = 0; i < D; ++i )
668   {
669     x_[ i ] /= a;
670   }
671   return *this;
672 }
673 
674 template < int D, class T >
675 inline bool Position< D, T >::operator==( const Position< D, T >& y ) const
676 {
677   for ( int i = 0; i < D; ++i )
678   {
679     if ( x_[ i ] != y.x_[ i ] )
680     {
681       return false;
682     }
683   }
684   return true;
685 }
686 
687 template < int D, class T >
688 inline bool Position< D, T >::operator!=( const Position< D, T >& y ) const
689 {
690   for ( int i = 0; i < D; ++i )
691   {
692     if ( x_[ i ] != y.x_[ i ] )
693     {
694       return true;
695     }
696   }
697   return false;
698 }
699 
700 template < int D, class T >
701 inline bool Position< D, T >::operator<( const Position< D, T >& y ) const
702 {
703   for ( int i = 0; i < D; ++i )
704   {
705     if ( x_[ i ] >= y.x_[ i ] )
706     {
707       return false;
708     }
709   }
710   return true;
711 }
712 
713 template < int D, class T >
714 inline bool Position< D, T >::operator>( const Position< D, T >& y ) const
715 {
716   for ( int i = 0; i < D; ++i )
717   {
718     if ( x_[ i ] <= y.x_[ i ] )
719     {
720       return false;
721     }
722   }
723   return true;
724 }
725 
726 template < int D, class T >
727 inline bool Position< D, T >::operator<=( const Position< D, T >& y ) const
728 {
729   for ( int i = 0; i < D; ++i )
730   {
731     if ( x_[ i ] > y.x_[ i ] )
732     {
733       return false;
734     }
735   }
736   return true;
737 }
738 
739 template < int D, class T >
740 inline bool Position< D, T >::operator>=( const Position< D, T >& y ) const
741 {
742   for ( int i = 0; i < D; ++i )
743   {
744     if ( x_[ i ] < y.x_[ i ] )
745     {
746       return false;
747     }
748   }
749   return true;
750 }
751 
752 template < int D, class T >
753 T
length()754 Position< D, T >::length() const
755 {
756   T lensq = 0;
757   for ( int i = 0; i < D; ++i )
758   {
759     lensq += x_[ i ] * x_[ i ];
760   }
761   return std::sqrt( lensq );
762 }
763 
764 template < int D, class T >
string()765 Position< D, T >::operator std::string() const
766 {
767   std::stringstream ss;
768   ss << *this;
769   return ss.str();
770 }
771 
772 template < int D, class T >
773 void
print(std::ostream & out,char sep)774 Position< D, T >::print( std::ostream& out, char sep ) const
775 {
776   out << x_[ 0 ];
777   for ( int i = 1; i < D; ++i )
778   {
779     out << sep << x_[ i ];
780   }
781 }
782 
783 template < int D, class T >
784 std::ostream& operator<<( std::ostream& os, const Position< D, T >& pos )
785 {
786   os << "(";
787   if ( D > 0 )
788   {
789     os << pos.x_[ 0 ];
790   }
791   for ( int i = 1; i < D; ++i )
792   {
793     os << ", " << pos.x_[ i ];
794   }
795   os << ")";
796   return os;
797 }
798 
799 } // namespace nest
800 
801 #endif
802