LORENE
des_prof_scalar.C
1/*
2 * Draws the profile of a {\tt Scalar} along some radial axis determined by
3 * a fixed value of $(\theta, \phi)$.
4 */
5
6/*
7 * Copyright (c) 2004-2005 Eric Gourgoulhon & Philippe Grandclement
8 *
9 * This file is part of LORENE.
10 *
11 * LORENE is free software; you can redistribute it and/or modify
12 * it under the terms of the GNU General Public License as published by
13 * the Free Software Foundation; either version 2 of the License, or
14 * (at your option) any later version.
15 *
16 * LORENE is distributed in the hope that it will be useful,
17 * but WITHOUT ANY WARRANTY; without even the implied warranty of
18 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
19 * GNU General Public License for more details.
20 *
21 * You should have received a copy of the GNU General Public License
22 * along with LORENE; if not, write to the Free Software
23 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
24 *
25 */
26
27
28char des_prof_scalar_C[] = "$Header: /cvsroot/Lorene/C++/Source/Non_class_members/Graphics/des_prof_scalar.C,v 1.12 2014/10/13 08:53:22 j_novak Exp $" ;
29
30
31/*
32 * $Id: des_prof_scalar.C,v 1.12 2014/10/13 08:53:22 j_novak Exp $
33 * $Log: des_prof_scalar.C,v $
34 * Revision 1.12 2014/10/13 08:53:22 j_novak
35 * Lorene classes and functions now belong to the namespace Lorene.
36 *
37 * Revision 1.11 2014/10/06 15:16:05 j_novak
38 * Modified #include directives to use c++ syntax.
39 *
40 * Revision 1.10 2012/01/17 10:35:40 j_penner
41 * added point plot
42 *
43 * Revision 1.9 2008/08/19 06:42:00 j_novak
44 * Minor modifications to avoid warnings with gcc 4.3. Most of them concern
45 * cast-type operations, and constant strings that must be defined as const char*
46 *
47 * Revision 1.8 2005/03/25 19:57:04 e_gourgoulhon
48 * Added plot of domain boundaries (new argument draw_bound).
49 *
50 * Revision 1.7 2004/05/17 19:47:25 e_gourgoulhon
51 * -- Function des_profile_mult(const Scalar**,...): added argument
52 * device.
53 * -- Functions des_meridian: added arguments device and closeit.
54 *
55 * Revision 1.6 2004/04/06 07:47:29 j_novak
56 * Added a #include "string.h"
57 *
58 * Revision 1.5 2004/04/05 14:42:02 e_gourgoulhon
59 * Added functions des_meridian.
60 *
61 * Revision 1.4 2004/02/17 22:18:00 e_gourgoulhon
62 * Changed prototype of des_profile_mult (added radial_scale, theta and
63 * phi can now vary from one profile to the other, etc...)
64 *
65 * Revision 1.3 2004/02/16 13:23:33 e_gourgoulhon
66 * Function des_profile_mult: added delete [] uutab at the end.
67 *
68 * Revision 1.2 2004/02/15 21:57:45 e_gourgoulhon
69 * des_profile_mult: changed argument Scalar* to Scalar**.
70 *
71 * Revision 1.1 2004/02/12 16:21:28 e_gourgoulhon
72 * Functions des_profile for Scalar's transfered from file des_prof_cmp.C.
73 * Added new function des_profile_mult.
74 *
75 *
76 * $Header: /cvsroot/Lorene/C++/Source/Non_class_members/Graphics/des_prof_scalar.C,v 1.12 2014/10/13 08:53:22 j_novak Exp $
77 *
78 */
79
80// Header C
81#include <cmath>
82#include <cstring>
83
84// Header Lorene
85#include "scalar.h"
86#include "graphique.h"
87
88#include <vector>
89namespace Lorene {
90//******************************************************************************
91// VERSION SCALAR SANS UNITES
92
93void des_profile(const Scalar& uu, double r_min, double r_max,
94 double theta, double phi, const char* nomy, const char* title,
95 bool draw_bound) {
96
97
98 const int npt = 400 ; // Number of points along the axis
99
100 float uutab[npt] ; // Value of uu at the npt points
101
102 double hr = (r_max - r_min) / double(npt-1) ;
103
104 for (int i=0; i<npt; i++) {
105
106 double r = hr * i + r_min ;
107
108 uutab[i] = float(uu.val_point(r, theta, phi)) ;
109
110 }
111
112 float xmin = float(r_min) ;
113 float xmax = float(r_max) ;
114
115 const char* nomx = "r" ;
116
117 if (title == 0x0) {
118 title = "" ;
119 }
120
121 if (nomy == 0x0) {
122 nomy = "" ;
123 }
124
125 // Preparations for the drawing of boundaries
126 // ------------------------------------------
127 const Map& mp = uu.get_mp() ;
128 int nz = mp.get_mg()->get_nzone() ;
129 int l_max = (mp.get_mg()->get_type_r(nz-1) == UNSURR) ? nz-2 : nz-1 ;
130
131 float* xbound = new float[l_max+1] ;
132 int nbound = 0 ;
133
134 if (draw_bound) {
135 const double xi_max = 1. ;
136 for (int l=0; l<=l_max; l++) {
137
138 double rb = mp.val_r(l, xi_max, theta, phi) ;
139
140 if ((rb >= r_min) && (rb <= r_max)) {
141 xbound[nbound] = float(rb) ;
142 nbound++ ;
143 }
144 }
145 }
146
147 des_profile(uutab, npt, xmin, xmax, nomx, nomy, title, 0x0,
148 nbound, xbound) ;
149
150 delete [] xbound ;
151
152}
153
154//******************************************************************************
155
156void des_profile(const Scalar& uu, double r_min, double r_max, double scale,
157 double theta, double phi, const char* nomx, const char* nomy,
158 const char* title, bool draw_bound) {
159
160
161 const int npt = 400 ; // Number of points along the axis
162
163 float uutab[npt] ; // Value of uu at the npt points
164
165 double hr = (r_max - r_min) / double(npt-1) ;
166
167 for (int i=0; i<npt; i++) {
168
169 double r = hr * i + r_min ;
170
171 uutab[i] = float(uu.val_point(r, theta, phi)) ;
172
173 }
174
175 float xmin = float(r_min * scale) ;
176 float xmax = float(r_max * scale) ;
177
178
179 if (title == 0x0) {
180 title = "" ;
181 }
182
183 if (nomx == 0x0) {
184 nomx = "" ;
185 }
186
187 if (nomy == 0x0) {
188 nomy = "" ;
189 }
190
191 // Preparations for the drawing of boundaries
192 // ------------------------------------------
193 const Map& mp = uu.get_mp() ;
194 int nz = mp.get_mg()->get_nzone() ;
195 int l_max = (mp.get_mg()->get_type_r(nz-1) == UNSURR) ? nz-2 : nz-1 ;
196
197 float* xbound = new float[l_max+1] ;
198 int nbound = 0 ;
199
200 if (draw_bound) {
201 const double xi_max = 1. ;
202 for (int l=0; l<=l_max; l++) {
203
204 double rb = mp.val_r(l, xi_max, theta, phi) ;
205
206 if ((rb >= r_min) && (rb <= r_max)) {
207 xbound[nbound] = float(rb) ;
208 nbound++ ;
209 }
210 }
211 }
212
213 // Call to the low level routine
214 // -----------------------------
215 des_profile(uutab, npt, xmin, xmax, nomx, nomy, title, 0x0,
216 nbound, xbound) ;
217
218 delete [] xbound ;
219
220}
221
222
223//******************************************************************************
224
225void des_profile_mult(const Scalar** uu, int nprof, double r_min, double r_max,
226 const double* theta, const double* phi, double radial_scale,
227 bool closeit, const char* nomy, const char* title, int ngraph,
228 const char* nomx, const int* line_style, const char* device,
229 bool draw_bound) {
230
231 // Special case of no graphical output:
232 if (device != 0x0) {
233 if ((device[0] == '/') && (device[1] == 'n')) return ;
234 }
235
236 const int npt = 400 ; // Number of points along the axis
237 double rr[npt] ;
238
239 float* uutab = new float[npt*nprof] ; // Value of uu at the npt points
240 // for each of the nprof profiles
241
242 double hr = (r_max - r_min) / double(npt-1) ;
243
244 for (int i=0; i<npt; i++) {
245 rr[i] = hr * i + r_min ;
246 }
247
248
249 for (int j=0; j<nprof; j++) {
250
251 const Scalar& vv = *(uu[j]) ;
252
253 for (int i=0; i<npt; i++) {
254 uutab[j*npt+i] = float(vv.val_point(rr[i], theta[j], phi[j])) ;
255 }
256 }
257
258
259 float xmin = float(radial_scale * r_min) ;
260 float xmax = float(radial_scale * r_max) ;
261
262 if (nomx == 0x0) nomx = "r" ;
263
264 if (nomy == 0x0) nomy = "" ;
265
266 if (title == 0x0) title = "" ;
267
268 // Preparations for the drawing of boundaries
269 // ------------------------------------------
270
271 int nbound_max = 100 * nprof ;
272 float* xbound = new float[nbound_max] ;
273 int nbound = 0 ;
274
275 if (draw_bound) {
276 const double xi_max = 1. ;
277 for (int j=0; j<nprof; j++) {
278
279 const Map& mp = uu[j]->get_mp() ;
280 int nz = mp.get_mg()->get_nzone() ;
281 int l_max = (mp.get_mg()->get_type_r(nz-1) == UNSURR) ? nz-2 : nz-1 ;
282
283 for (int l=0; l<=l_max; l++) {
284
285 double rb = mp.val_r(l, xi_max, theta[j], phi[j]) ;
286
287 if ((rb >= r_min) && (rb <= r_max)) {
288 xbound[nbound] = float(rb * radial_scale) ;
289 nbound++ ;
290 if (nbound > nbound_max-1) {
291 cout << "des_profile_mult : nbound too large !" << endl ;
292 abort() ;
293 }
294 }
295 }
296 }
297 }
298
299 // Call to the low level routine
300 // -----------------------------
301
302 des_profile_mult(uutab, nprof, npt, xmin, xmax, nomx, nomy, title,
303 line_style, ngraph, closeit, device, nbound, xbound) ;
304
305
306 delete [] uutab ;
307 delete [] xbound ;
308
309}
310
311//******************************************************************************
312
313void des_meridian(const Scalar& uu, double r_min, double r_max,
314 const char* nomy, int ngraph, const char* device,
315 bool closeit, bool draw_bound) {
316
317 // Special case of no graphical output:
318 if (device != 0x0) {
319 if ((device[0] == '/') && (device[1] == 'n')) return ;
320 }
321
322 const Scalar* des[] = {&uu, &uu, &uu, &uu, &uu} ;
323 double phi1[] = {0., 0., 0., 0.25*M_PI, 0.25*M_PI} ;
324 double theta1[] = {0., 0.25*M_PI, 0.5*M_PI, 0., 0.25*M_PI} ;
325
326 des_profile_mult(des, 5, r_min, r_max, theta1, phi1, 1., closeit,
327 nomy,
328 "phi=0: th=0, pi/4, pi/2, phi=pi/4: th=0, pi/4",
329 ngraph, 0x0, 0x0, device, draw_bound) ;
330
331}
332
333//******************************************************************************
334
335
336void des_meridian(const Sym_tensor& hh, double r_min, double r_max,
337 const char* name, int ngraph0, const char* device,
338 bool closeit) {
339
340 // Special case of no graphical output:
341 if (device != 0x0) {
342 if ((device[0] == '/') && (device[1] == 'n')) return ;
343 }
344
345 char nomy[80] ;
346
347 int k = 0 ;
348 for (int i=1; i<=3; i++) {
349 for (int j=i; j<=3; j++) {
350
351 char nom_i[3] ;
352 sprintf(nom_i, "%d", i) ;
353 char nom_j[3] ;
354 sprintf(nom_j, "%d", j) ;
355 strncpy(nomy, name, 40) ;
356 strcat(nomy, " comp. ") ;
357 strcat(nomy, nom_i) ;
358 strcat(nomy, nom_j) ;
359
360 des_meridian(hh(i,j), r_min, r_max, nomy, ngraph0+k, device,
361 closeit) ;
362 k++ ;
363
364 }
365 }
366
367}
368
369
370//******************************************************************************
371// VERSION SCALAR SANS UNITES
372
373void des_points(const Scalar& uu,
374 double theta, double phi, const char* nomy, const char* title,
375 bool draw_bound) {
376
377 const Map& mp = uu.get_mp() ;
378 int nz = mp.get_mg()->get_nzone() ;
379 int nt = mp.get_mg()->get_nt(nz-1) ;
380 int np = mp.get_mg()->get_np(nz-1) ;
381
382// const int npt = *(uu.get_mp().get_mg())->get_nzone() ; // Number of points along the axis
383
384
385 int npt=0;
386
387 for(int ii = 0; ii<nz; ii++)
388 npt += (uu.get_mp().get_mg())->get_nr(ii) ;
389
390 float *uutab = new float[npt] ; // define a dynamic array
391 float *xtab = new float[npt] ; // define a dynamic array
392
393 Mtbl r = *(mp.r.c);
394
395 for(int ii = 0; ii<nz; ii++){
396 int nr = (uu.get_mp().get_mg())->get_nr(ii) ;
397 for(int ij=0; ij<nr; ij++){
398 uutab[ii*nr+ij] = float(uu.val_grid_point(ii,np-1,nt-1,ij)) ;
399 xtab[ii*nr+ij] = float(r(ii,np-1,nt-1,ij)) ;
400 }
401 }
402
403 float xmin = float(totalmin(r)) ;
404 float xmax = float(totalmax(r)) ;
405
406 const char* nomx = "r" ;
407
408 if (title == 0x0) {
409 title = "" ;
410 }
411
412 if (nomy == 0x0) {
413 nomy = "" ;
414 }
415
416 // Preparations for the drawing of boundaries
417 // ------------------------------------------
418 int l_max = (mp.get_mg()->get_type_r(nz-1) == UNSURR) ? nz-2 : nz-1 ;
419
420 float* xbound = new float[l_max+1] ;
421 int nbound = 0 ;
422
423 if (draw_bound) {
424 const double xi_max = 1. ;
425 for (int l=0; l<=l_max; l++) {
426
427 double rb = mp.val_r(l, xi_max, theta, phi) ;
428
429 if ((rb >= xmin) && (rb <= xmax)) {
430 xbound[nbound] = float(rb) ;
431 nbound++ ;
432 }
433 }
434 }
435
436 des_profile(uutab, npt, xtab, nomx, nomy, title, 0x0,
437 nbound, xbound) ;
438
439 delete [] xbound ;
440
441}
442
443//******************************************************************************
444
445void des_points(const Scalar& uu, double scale,
446 double theta, double phi, const char* nomx, const char* nomy,
447 const char* title, bool draw_bound) {
448
449 const Map& mp = uu.get_mp() ;
450 int nz = mp.get_mg()->get_nzone() ;
451 int nt = mp.get_mg()->get_nt(nz-1) ;
452 int np = mp.get_mg()->get_np(nz-1) ;
453
454// const int npt = *(uu.get_mp().get_mg())->get_nzone() ; // Number of points along the axis
455
456
457 int npt=0;
458
459 for(int ii = 0; ii<nz; ii++)
460 npt += (uu.get_mp().get_mg())->get_nr(ii) ;
461
462 float *uutab = new float[npt] ; // define a dynamic array
463 float *xtab = new float[npt] ; // define a dynamic array
464
465 Mtbl r = *(mp.r.c);
466
467 for(int ii = 0; ii<nz; ii++){
468 int nr = (uu.get_mp().get_mg())->get_nr(ii) ;
469 for(int ij=0; ij<nr; ij++){
470 uutab[ii*nr+ij] = float(uu.val_grid_point(ii,np-1,nt-1,ij)) ;
471 xtab[ii*nr+ij] = float(r(ii,np-1,nt-1,ij)) ;
472 }
473 }
474
475 float xmin = float(totalmin(r) * scale) ;
476 float xmax = float(totalmax(r) * scale) ;
477
478
479 if (title == 0x0) {
480 title = "" ;
481 }
482
483 if (nomx == 0x0) {
484 nomx = "" ;
485 }
486
487 if (nomy == 0x0) {
488 nomy = "" ;
489 }
490
491 // Preparations for the drawing of boundaries
492 // ------------------------------------------
493 int l_max = (mp.get_mg()->get_type_r(nz-1) == UNSURR) ? nz-2 : nz-1 ;
494
495 float* xbound = new float[l_max+1] ;
496 int nbound = 0 ;
497
498 if (draw_bound) {
499 const double xi_max = 1. ;
500 for (int l=0; l<=l_max; l++) {
501
502 double rb = mp.val_r(l, xi_max, theta, phi) ;
503
504 if ((rb >= xmin/scale) && (rb <= xmax/scale)) {
505 xbound[nbound] = float(rb) ;
506 nbound++ ;
507 }
508 }
509 }
510
511 // Call to the low level routine
512 // -----------------------------
513 des_profile(uutab, npt, xtab, nomx, nomy, title, 0x0,
514 nbound, xbound) ;
515
516 delete [] xbound ;
517
518}
519}
void des_profile_mult(const float *uutab, int nprof, int nx, float xmin, float xmax, const char *nomx, const char *nomy, const char *title, const int *line_style, int ngraph, bool closeit, const char *device=0x0, int nbound=0, float *xbound=0x0)
Basic routine for drawing multiple profiles with uniform x sampling.
void des_profile(const float *uutab, int nx, float xmin, float xmax, const char *nomx, const char *nomy, const char *title, const char *device=0x0, int nbound=0, float *xbound=0x0)
Basic routine for drawing a single profile with uniform x sampling.
void des_points(const float *uutab, int nx, float xmin, float xmax, const char *nomx=0x0, const char *nomy=0x0, const char *title=0x0, const char *device=0x0, int nbound=0, float *xbound=0x0)
Basic routine for plotting points using grid locations.
void des_meridian(const Scalar &uu, double r_min, double r_max, const char *nomy, int ngraph, const char *device=0x0, bool closeit=false, bool draw_bound=true)
Draws 5 profiles of a scalar field along various radial axes in two meridional planes and .
double totalmin(const Mtbl &)
Minimum value of the Mtbl elements in all domain.
Definition mtbl_math.C:522
double totalmax(const Mtbl &)
Maximum value of the Mtbl elements in all domains.
Definition mtbl_math.C:494
Lorene prototypes.
Definition app_hor.h:64