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