LORENE
map.C
1 /*
2  * Methods of class Map
3  *
4  * (see file map.h for documentation)
5  *
6  */
7 
8 /*
9  * Copyright (c) 1999-2003 Eric Gourgoulhon
10  *
11  * This file is part of LORENE.
12  *
13  * LORENE is free software; you can redistribute it and/or modify
14  * it under the terms of the GNU General Public License as published by
15  * the Free Software Foundation; either version 2 of the License, or
16  * (at your option) any later version.
17  *
18  * LORENE is distributed in the hope that it will be useful,
19  * but WITHOUT ANY WARRANTY; without even the implied warranty of
20  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
21  * GNU General Public License for more details.
22  *
23  * You should have received a copy of the GNU General Public License
24  * along with LORENE; if not, write to the Free Software
25  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
26  *
27  */
28 
29 
30 char map_C[] = "$Header: /cvsroot/Lorene/C++/Source/Map/map.C,v 1.11 2014/10/13 08:53:02 j_novak Exp $" ;
31 
32 
33 /*
34  * $Id: map.C,v 1.11 2014/10/13 08:53:02 j_novak Exp $
35  * $Log: map.C,v $
36  * Revision 1.11 2014/10/13 08:53:02 j_novak
37  * Lorene classes and functions now belong to the namespace Lorene.
38  *
39  * Revision 1.10 2014/10/06 15:13:12 j_novak
40  * Modified #include directives to use c++ syntax.
41  *
42  * Revision 1.9 2008/09/29 13:23:51 j_novak
43  * Implementation of the angular mapping associated with an affine
44  * mapping. Things must be improved to take into account the domain index.
45  *
46  * Revision 1.8 2004/01/29 08:50:03 p_grandclement
47  * Modification of Map::operator==(const Map&) and addition of the surface
48  * integrales using Scalar.
49  *
50  * Revision 1.7 2003/12/30 22:53:23 e_gourgoulhon
51  * Added methods flat_met_spher() and flat_met_cart() to get
52  * flat metric associated with the coordinates described by the mapping.
53  *
54  * Revision 1.6 2003/10/15 10:30:46 e_gourgoulhon
55  * Map::set_ori: changed x,y,z to xo,yo,zo not to shadow Coord's x,y,z
56  * Map::set_rot_phi: changed phi to newphi not to shadow Coord phi.
57  *
58  * Revision 1.5 2003/06/20 14:45:06 f_limousin
59  * Add operator==
60  *
61  * Revision 1.4 2002/10/16 14:36:40 j_novak
62  * Reorganization of #include instructions of standard C++, in order to
63  * use experimental version 3 of gcc.
64  *
65  * Revision 1.3 2002/05/22 12:44:04 f_limousin
66  * Added print of ori_x, ori_y, ori_z and rot_phi in operator<<.
67  *
68  * Revision 1.2 2001/12/04 21:27:53 e_gourgoulhon
69  *
70  * All writing/reading to a binary file are now performed according to
71  * the big endian convention, whatever the system is big endian or
72  * small endian, thanks to the functions fwrite_be and fread_be
73  *
74  * Revision 1.1.1.1 2001/11/20 15:19:27 e_gourgoulhon
75  * LORENE
76  *
77  * Revision 2.12 2000/02/11 14:26:54 eric
78  * *** empty log message ***
79  *
80  * Revision 2.11 2000/02/11 13:38:23 eric
81  * Ajout de la fonction convert_absolute.
82  *
83  * Revision 2.10 2000/02/09 13:25:29 eric
84  * Mise en conformite avec le nouveau constructeur de Base_vect_spher
85  * (passage des arguments ori_x, ori_y et ori_z).
86  *
87  * Revision 2.9 2000/01/12 12:54:53 eric
88  * Ajout du Cmp null, *p_cmp_zero, et de la methode associee cmp_zero().
89  *
90  * Revision 2.8 2000/01/10 13:28:47 eric
91  * Ajout des bases vectorielles associees aux coordonnees :
92  * membres bvect_spher et bvect_cart.
93  *
94  * Revision 2.7 1999/11/24 14:39:34 eric
95  * Appel de reset_coord() dans les fonctions set_ori et set_rot_phi.
96  *
97  * Revision 2.6 1999/11/22 10:35:36 eric
98  * Fonction del_coord() rebaptisee reset_coord().
99  *
100  * Revision 2.5 1999/10/15 14:12:29 eric
101  * Suppression de l'appel a del_coord() dans le destructeur de Map.
102  *
103  * Revision 2.4 1999/10/15 09:23:04 eric
104  * *** empty log message ***
105  *
106  * Revision 2.3 1999/10/14 14:26:51 eric
107  * Depoussierage.
108  * Documentation.
109  *
110  * Revision 2.2 1999/03/01 16:57:13 eric
111  * Operateur <<
112  *
113  * Revision 2.1 1999/03/01 14:57:52 eric
114  * *** empty log message ***
115  *
116  * Revision 2.0 1999/02/15 10:42:45 hyc
117  * *** empty log message ***
118  *
119  * $Header: /cvsroot/Lorene/C++/Source/Map/map.C,v 1.11 2014/10/13 08:53:02 j_novak Exp $
120  *
121  */
122 
123 // headers C
124 #include <cmath>
125 
126 // headers Lorene
127 #include "map.h"
128 #include "cmp.h"
129 #include "metric.h"
130 #include "utilitaires.h"
131 
132  //---------------//
133  // Constructeurs //
134  //---------------//
135 
136 // Constructor from a grid
137 // -----------------------
138 namespace Lorene {
139 Map::Map(const Mg3d& mgi) : mg(&mgi),
140  ori_x(0), ori_y(0), ori_z(0), rot_phi(0),
141  bvect_spher(ori_x, ori_y, ori_z, rot_phi,
142  "Mapping orthonormal spherical basis"),
143  bvect_cart(rot_phi, "Mapping Cartesian basis"),
144  p_flat_met_spher(0x0),
145  p_flat_met_cart(0x0),
146  p_mp_angu(0x0)
147 {
148  // The Coord's are constructed by the default Coord constructor
149 
150  // The null Cmp :
151  p_cmp_zero = new Cmp(this) ;
153 }
154 
155 // Copy constructor
156 // ----------------
157 Map::Map(const Map& mp) : mg(mp.mg),
158  ori_x(mp.ori_x), ori_y(mp.ori_y), ori_z(mp.ori_z),
159  rot_phi(mp.rot_phi),
160  bvect_spher(ori_x, ori_y, ori_z, rot_phi,
161  "Mapping orthonormal spherical basis"),
162  bvect_cart(rot_phi, "Mapping Cartesian basis"),
163  p_flat_met_spher(0x0),
164  p_flat_met_cart(0x0),
165  p_mp_angu(0x0)
166 {
167  // The Coord's are constructed by the default Coord constructor
168 
169  // The null Cmp :
170  p_cmp_zero = new Cmp(this) ;
172 }
173 
174 // Constructor from file
175 // ---------------------
176 Map::Map(const Mg3d& mgi, FILE* fd) : mg(&mgi),
177  bvect_spher(0., 0., 0., 0.,
178  "Mapping orthonormal spherical basis"),
179  bvect_cart(0., "Mapping Cartesian basis"),
180  p_flat_met_spher(0x0),
181  p_flat_met_cart(0x0),
182  p_mp_angu(0x0)
183 {
184  Mg3d* mg_tmp = new Mg3d(fd) ; // la multi-grille d'origine
185  if (*mg != *mg_tmp) {
186  cout << "Map::Map(const Mg3d&, FILE*): grid not consistent !"
187  << endl ;
188  abort() ;
189  }
190  delete mg_tmp ;
191 
192  fread_be(&ori_x, sizeof(double), 1, fd) ; // ori_x
193  fread_be(&ori_y, sizeof(double), 1, fd) ; // ori_y
194  fread_be(&ori_z, sizeof(double), 1, fd) ; // ori_z
195  fread_be(&rot_phi, sizeof(double), 1, fd) ; // rot_phi
196 
200 
201  // The Coord's are constructed by the default Coord constructor
202 
203  // The null Cmp :
204  p_cmp_zero = new Cmp(this) ;
206 }
207 
208  //--------------//
209  // Destructeurs //
210  //--------------//
211 
212 // Destructeur
214  if (p_flat_met_spher != 0x0) delete p_flat_met_spher ;
215  if (p_flat_met_cart != 0x0) delete p_flat_met_cart ;
216  if (p_mp_angu != 0x0) delete p_mp_angu ;
217  delete p_cmp_zero ;
218 }
219 
220  //------------//
221  // Sauvegarde //
222  //------------//
223 
224 void Map::sauve(FILE* fd) const {
225 
226  mg->sauve(fd) ; // la multi-grille
227 
228  fwrite_be(&ori_x, sizeof(double), 1, fd) ; // ori_x
229  fwrite_be(&ori_y, sizeof(double), 1, fd) ; // ori_y
230  fwrite_be(&ori_z, sizeof(double), 1, fd) ; // ori_z
231  fwrite_be(&rot_phi, sizeof(double), 1, fd) ; // rot_phi
232 }
233 
234  //------------//
235  // Impression //
236  //------------//
237 
238 // Operateurs <<
239 ostream& operator<<(ostream& o, const Map & cv) {
240  o << "Absolute coordinates of the mapping origin: " << endl ;
241  o << " X_0, Y_0, Z_0 : " << cv.get_ori_x() << " "
242  << cv.get_ori_y() << " " << cv.get_ori_z() << endl ;
243  o << "Rotation angle between the x-axis and X-axis : "
244  << cv.get_rot_phi() << endl ;
245  cv >> o ;
246  return o ;
247 }
248 
249  //-------------------//
250  // Methodes diverses //
251  //-------------------//
252 
253 void Map::set_ori(double xo, double yo, double zo) {
254  ori_x = xo ;
255  ori_y = yo ;
256  ori_z = zo ;
257 
259 
260  reset_coord() ; // Mise a jour des Coords
261 }
262 
263 void Map::set_rot_phi(double newphi) {
264 
265  rot_phi = newphi ;
266 
269 
270  reset_coord() ; // Mise a jour des Coords
271 
272 }
273 
274 // Gestion de la memoire
275 // ---------------------
277 
278  // Coordonnees spheriques centrees sur la grille :
279  r.del_t() ;
280  tet.del_t() ;
281  phi.del_t() ;
282  sint.del_t() ;
283  cost.del_t() ;
284  sinp.del_t() ;
285  cosp.del_t() ;
286 
287  // Coordonnees cartesiennes centrees sur la grille :
288  x.del_t() ;
289  y.del_t() ;
290  z.del_t() ;
291 
292  // Coordonnees cartesiennes "absolues" :
293  xa.del_t() ;
294  ya.del_t() ;
295  za.del_t() ;
296 }
297 
298 
299 // Conversion coordonnees (X,Y,Z) --> (r, theta, phi)
300 // --------------------------------------------------
301 
302 void Map::convert_absolute(double xx, double yy, double zz,
303  double& rr, double& theta, double& pphi) const {
304 
305  // Cartesian coordinates aligned with the absolute ones but centered
306  // on the mapping origin :
307  double x1 = xx - ori_x ;
308  double y1 = yy - ori_y ;
309  double z1 = zz - ori_z ;
310 
311  // Spherical coordinates :
312  double rho2 = x1*x1 + y1*y1 ;
313  double rho = sqrt( rho2 ) ;
314  rr = sqrt(rho2 + z1*z1) ;
315  theta = atan2(rho, z1) ;
316  pphi = atan2(y1, x1) - rot_phi ; // (rotation)
317  if (pphi < 0) pphi += 2*M_PI ;
318 
319 }
320 
322 
323  if (p_flat_met_spher == 0x0) {
324  p_flat_met_spher = new Metric_flat(*this, bvect_spher) ;
325  }
326 
327  return *p_flat_met_spher ;
328 
329 }
330 
332 
333  if (p_flat_met_cart == 0x0) {
334  p_flat_met_cart = new Metric_flat(*this, bvect_cart) ;
335  }
336 
337  return *p_flat_met_cart ;
338 
339 }
340 
341 
342 
343 
344 
345 
346 
347 
348 
349 
350 
351 }
Lorene::Map::z
Coord z
z coordinate centered on the grid
Definition: map.h:728
Lorene::Map::set_rot_phi
void set_rot_phi(double phi0)
Sets a new rotation angle.
Definition: map.C:263
Lorene::Map::p_cmp_zero
Cmp * p_cmp_zero
The null Cmp.
Definition: map.h:713
Lorene::Map::get_ori_x
double get_ori_x() const
Returns the x coordinate of the origin.
Definition: map.h:768
Lorene::Map::cost
Coord cost
Definition: map.h:722
Lorene::Mg3d
Multi-domain grid.
Definition: grilles.h:273
Lorene::Base_vect_spher::set_rot_phi
void set_rot_phi(double rot_phi_i)
Sets a new value to the angle rot_phi between the x –axis and the absolute frame X –axis.
Definition: base_vect_spher.C:184
Lorene::Map::p_mp_angu
Map_af * p_mp_angu
Pointer on the "angular" mapping.
Definition: map.h:715
Lorene::Map::get_rot_phi
double get_rot_phi() const
Returns the angle between the x –axis and X –axis.
Definition: map.h:775
Lorene
Lorene prototypes.
Definition: app_hor.h:64
Lorene::Map::sauve
virtual void sauve(FILE *) const
Save in a file.
Definition: map.C:224
Lorene::Map::set_ori
void set_ori(double xa0, double ya0, double za0)
Sets a new origin.
Definition: map.C:253
Lorene::Map::flat_met_spher
const Metric_flat & flat_met_spher() const
Returns the flat metric associated with the spherical coordinates and with components expressed in th...
Definition: map.C:321
Lorene::Map::get_ori_y
double get_ori_y() const
Returns the y coordinate of the origin.
Definition: map.h:770
Lorene::Base_vect_cart::set_rot_phi
void set_rot_phi(double rot_phi_i)
Sets a new value to the angle rot_phi between the x –axis and the absolute frame X –axis.
Definition: base_vect_cart.C:212
Lorene::Map::sint
Coord sint
Definition: map.h:721
Lorene::Map::r
Coord r
r coordinate centered on the grid
Definition: map.h:718
Lorene::Map::p_flat_met_spher
Metric_flat * p_flat_met_spher
Pointer onto the flat metric associated with the spherical coordinates and with components expressed ...
Definition: map.h:702
Lorene::Cmp::set_etat_zero
void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Definition: cmp.C:289
Lorene::Metric_flat
Flat metric for tensor calculation.
Definition: metric.h:261
Lorene::Map::sinp
Coord sinp
Definition: map.h:723
Lorene::Map::p_flat_met_cart
Metric_flat * p_flat_met_cart
Pointer onto the flat metric associated with the Cartesian coordinates and with components expressed ...
Definition: map.h:707
Lorene::Cmp
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Definition: cmp.h:446
Lorene::Map::ori_y
double ori_y
Absolute coordinate y of the origin.
Definition: map.h:679
Lorene::Mg3d::sauve
void sauve(FILE *fd, bool save_base=false) const
Saves into a file.
Definition: mg3d.C:371
Lorene::Map::flat_met_cart
const Metric_flat & flat_met_cart() const
Returns the flat metric associated with the Cartesian coordinates and with components expressed in th...
Definition: map.C:331
Lorene::Map::ori_x
double ori_x
Absolute coordinate x of the origin.
Definition: map.h:678
Lorene::Map::tet
Coord tet
coordinate centered on the grid
Definition: map.h:719
Lorene::fread_be
int fread_be(int *aa, int size, int nb, FILE *fich)
Reads integer(s) from a binary file according to the big endian convention.
Definition: fread_be.C:69
Lorene::Map::phi
Coord phi
coordinate centered on the grid
Definition: map.h:720
Lorene::sqrt
Cmp sqrt(const Cmp &)
Square root.
Definition: cmp_math.C:220
Lorene::Map::bvect_cart
Base_vect_cart bvect_cart
Cartesian basis associated with the coordinates (x,y,z) of the mapping, i.e.
Definition: map.h:697
Lorene::Map::x
Coord x
x coordinate centered on the grid
Definition: map.h:726
Lorene::Map::xa
Coord xa
Absolute x coordinate.
Definition: map.h:730
Lorene::Map::reset_coord
virtual void reset_coord()
Resets all the member Coords.
Definition: map.C:276
Lorene::Base_vect_spher::set_ori
void set_ori(double xa0, double ya0, double za0)
Sets a new origin.
Definition: base_vect_spher.C:175
Lorene::Map::bvect_spher
Base_vect_spher bvect_spher
Orthonormal vectorial basis associated with the coordinates of the mapping.
Definition: map.h:689
Lorene::Map::get_ori_z
double get_ori_z() const
Returns the z coordinate of the origin.
Definition: map.h:772
Lorene::Map::~Map
virtual ~Map()
Destructor.
Definition: map.C:213
Lorene::Map::ori_z
double ori_z
Absolute coordinate z of the origin.
Definition: map.h:680
Lorene::Map::convert_absolute
void convert_absolute(double xx, double yy, double zz, double &rr, double &theta, double &pphi) const
Determines the coordinates corresponding to given absolute Cartesian coordinates (X,...
Definition: map.C:302
Lorene::fwrite_be
int fwrite_be(const int *aa, int size, int nb, FILE *fich)
Writes integer(s) into a binary file according to the big endian convention.
Definition: fwrite_be.C:70
Lorene::Map::rot_phi
double rot_phi
Angle between the x –axis and X –axis.
Definition: map.h:681
Lorene::Coord::del_t
void del_t() const
Logical destructor (deletes the Mtbl member *c ).
Definition: coord.C:125
Lorene::Map::y
Coord y
y coordinate centered on the grid
Definition: map.h:727
Lorene::Map::za
Coord za
Absolute z coordinate.
Definition: map.h:732
Lorene::Map::cosp
Coord cosp
Definition: map.h:724
Lorene::Map
Base class for coordinate mappings.
Definition: map.h:670
Lorene::Map::ya
Coord ya
Absolute y coordinate.
Definition: map.h:731
Lorene::Map::Map
Map(const Mg3d &)
Constructor from a multi-domain 3D grid.
Definition: map.C:139
Lorene::Map::mg
const Mg3d * mg
Pointer on the multi-grid Mgd3 on which this is defined
Definition: map.h:676