LORENE
scalar_visu.C
1/*
2 * 3D visualization of a Scalar via OpenDX
3 *
4 * (see file scalar.h for documentation).
5 *
6 */
7
8/*
9 * Copyright (c) 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 version 2
15 * as published by the Free Software Foundation.
16 *
17 * LORENE is distributed in the hope that it will be useful,
18 * but WITHOUT ANY WARRANTY; without even the implied warranty of
19 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
20 * GNU General Public License for more details.
21 *
22 * You should have received a copy of the GNU General Public License
23 * along with LORENE; if not, write to the Free Software
24 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
25 *
26 */
27
28
29
30/*
31 * $Id: scalar_visu.C,v 1.10 2016/12/05 16:18:19 j_novak Exp $
32 * $Log: scalar_visu.C,v $
33 * Revision 1.10 2016/12/05 16:18:19 j_novak
34 * Suppression of some global variables (file names, loch, ...) to prevent redefinitions
35 *
36 * Revision 1.9 2014/10/13 08:53:47 j_novak
37 * Lorene classes and functions now belong to the namespace Lorene.
38 *
39 * Revision 1.8 2014/10/06 15:16:16 j_novak
40 * Modified #include directives to use c++ syntax.
41 *
42 * Revision 1.7 2005/02/18 13:14:19 j_novak
43 * Changing of malloc/free to new/delete + suppression of some unused variables
44 * (trying to avoid compilation warnings).
45 *
46 * Revision 1.6 2004/03/11 12:07:55 e_gourgoulhon
47 * Added method visu_section_anim.
48 *
49 * Revision 1.5 2003/12/19 15:18:17 j_novak
50 * Shadow variables hunt
51 *
52 * Revision 1.4 2003/12/16 06:32:57 e_gourgoulhon
53 * Added method visu_box.
54 *
55 * Revision 1.3 2003/12/15 08:30:40 p_grandclement
56 * Addition of #include <string.h>
57 *
58 * Revision 1.2 2003/12/14 21:49:14 e_gourgoulhon
59 * Added argument start_dx (to launch OpenDX as a subprocess).
60 *
61 * Revision 1.1 2003/12/11 16:20:25 e_gourgoulhon
62 * First version.
63 *
64 *
65 * $Header: /cvsroot/Lorene/C++/Source/Tensor/Scalar/scalar_visu.C,v 1.10 2016/12/05 16:18:19 j_novak Exp $
66 *
67 */
68
69// C headers
70#include <cstdlib>
71#include <cstring>
72
73// Lorene headers
74#include "tensor.h"
75
76 //-----------------------------------------//
77 // visu_section : special cases //
78 //-----------------------------------------//
79
80namespace Lorene {
81void Scalar::visu_section(const char section_type, double aa, double umin,
82 double umax, double vmin, double vmax, const char* title0,
83 const char* filename0, bool start_dx, int nu, int nv) const {
84
85 Tbl plane(3,3) ;
86 plane.set_etat_qcq() ; // Memory allocation for the Tbl
87
88 switch (section_type) {
89
90 case 'x' : {
91 plane.set(0,0) = aa ; // Origin in the plane
92 plane.set(0,1) = 0. ; // (absolute Cartesian coordinates)
93 plane.set(0,2) = 0. ; //
94
95 plane.set(1,0) = 0. ; // u-coordinate unit vector
96 plane.set(1,1) = 1. ; // (absolute Cartesian components)
97 plane.set(1,2) = 0. ;
98
99 plane.set(2,0) = 0. ; // v-coordinate unit vector
100 plane.set(2,1) = 0. ; // (absolute Cartesian components)
101 plane.set(2,2) = 1. ;
102
103 visu_section(plane, umin, umax, vmin, vmax, title0, filename0,
104 start_dx, nu, nv) ;
105 break ;
106 }
107
108 case 'y' : {
109 plane.set(0,0) = 0. ; // Origin in the plane
110 plane.set(0,1) = aa ; // (absolute Cartesian coordinates)
111 plane.set(0,2) = 0. ; //
112
113 plane.set(1,0) = 1. ; // u-coordinate unit vector
114 plane.set(1,1) = 0. ; // (absolute Cartesian components)
115 plane.set(1,2) = 0. ;
116
117 plane.set(2,0) = 0. ; // v-coordinate unit vector
118 plane.set(2,1) = 0. ; // (absolute Cartesian components)
119 plane.set(2,2) = 1. ;
120
121 visu_section(plane, umin, umax, vmin, vmax, title0, filename0,
122 start_dx, nu, nv) ;
123 break ;
124 }
125
126 case 'z' : {
127 plane.set(0,0) = 0. ; // Origin in the plane
128 plane.set(0,1) = 0. ; // (absolute Cartesian coordinates)
129 plane.set(0,2) = aa ; //
130
131 plane.set(1,0) = 1. ; // u-coordinate unit vector
132 plane.set(1,1) = 0. ; // (absolute Cartesian components)
133 plane.set(1,2) = 0. ;
134
135 plane.set(2,0) = 0. ; // v-coordinate unit vector
136 plane.set(2,1) = 1. ; // (absolute Cartesian components)
137 plane.set(2,2) = 0. ;
138
139 visu_section(plane, umin, umax, vmin, vmax, title0, filename0,
140 start_dx, nu, nv) ;
141 break ;
142 }
143
144 default : {
145 cerr << "Scalar::visu_section : unknown type of section ! \n" ;
146 cerr << " section_type = " << section_type << endl ;
147 break ;
148 }
149 }
150
151}
152
153
154 //-----------------------------------------//
155 // visu_section : general case //
156 //-----------------------------------------//
157
158
159void Scalar::visu_section(const Tbl& plane, double umin, double umax,
160 double vmin, double vmax, const char* title0, const char* filename0,
161 bool start_dx, int nu, int nv) const {
162
163 char* title ;
164 char* title_quotes ;
165 if (title0 == 0x0) {
166 title = new char[2] ;
167 strcpy(title, " ") ;
168
169 title_quotes = new char[4] ;
170 strcpy(title_quotes, "\" \"") ;
171 }
172 else {
173 title = new char[ strlen(title0)+1 ] ;
174 strcpy(title, title0) ;
175
176 title_quotes = new char[ strlen(title0)+3 ] ;
177 strcpy(title_quotes, "\"") ;
178 strcat(title_quotes, title0) ;
179 strcat(title_quotes, "\"") ;
180 }
181
182 // --------------------------------------------------------
183 // Data file for OpenDX
184 // --------------------------------------------------------
185
186 char* filename ;
187 if (filename0 == 0x0) {
188 filename = new char[30] ;
189 strcpy(filename, "scalar_section.dxdata") ;
190 }
191 else {
192 filename = new char[ strlen(filename0)+8 ] ;
193 strcpy(filename, filename0) ;
194 strcat(filename, ".dxdata") ;
195 }
196
197 ofstream fdata(filename) ; // output file
198
199 fdata << title << "\n" ;
200 fdata << "size : " << nu << " x " << nv << "\n" ;
201 fdata << "u_min = " << umin << " u_max = " << umax << "\n" ;
202 fdata << "v_min = " << vmin << " v_max = " << vmax << "\n" ;
203
204 // Plane characterization
205 // ----------------------
206
207 double xa0 = plane(0,0) ;
208 double ya0 = plane(0,1) ;
209 double za0 = plane(0,2) ;
210
211 double eux = plane(1,0) ;
212 double euy = plane(1,1) ;
213 double euz = plane(1,2) ;
214
215 double evx = plane(2,0) ;
216 double evy = plane(2,1) ;
217 double evz = plane(2,2) ;
218
219
220 // The spectral coefficients are required
221 va.coef() ;
222 const Mtbl_cf& cva = *(va.c_cf) ;
223
224 // What follows assumes that the mapping is radial:
225 assert( dynamic_cast<const Map_radial*>(mp) != 0x0 ) ;
226
227 fdata.precision(5) ;
228 fdata.setf(ios::scientific) ;
229
230 // Loop on the points in the section plane
231 // ---------------------------------------
232 double du = (umax - umin) / double(nu-1) ;
233 double dv = (vmax - vmin) / double(nv-1) ;
234
235 int npoint = 0 ; // number of data points per line in the file
236
237 for (int j=0; j<nv; j++) {
238
239 double v = vmin + dv * j ;
240
241 for (int i=0; i<nu; i++) {
242
243 double u = umin + du * i ;
244
245 double xa = xa0 + u * eux + v * evx ;
246 double ya = ya0 + u * euy + v * evy ;
247 double za = za0 + u * euz + v * evz ;
248
249
250 // Values of (r,theta,phi) corresponding to (xa,ya,za) :
251 double rr, th, ph ; // polar coordinates of the mapping associated
252 // to *this
253
254 mp->convert_absolute(xa, ya, za, rr, th, ph) ;
255
256 // Values of (l,xi,theta',phi') corresponding to (r,theta,phi):
257 double xi ;
258 int l ;
259
260 mp->val_lx(rr, th, ph, l, xi) ; // radial mapping assumption
261
262 // Field value at this point:
263
264 double ff = cva.val_point(l, xi, th, ph) ;
265
266 fdata.width(14) ;
267 fdata << ff ;
268 npoint++ ;
269
270 if (npoint == 6) {
271 fdata << "\n" ;
272 npoint = 0 ;
273 }
274
275 }
276 }
277
278 if (npoint != 0) fdata << "\n" ;
279
280 fdata.close() ;
281
282 // --------------------------------------------------------
283 // Header file for OpenDX
284 // --------------------------------------------------------
285
286 char* headername ;
287 if (filename0 == 0x0) {
288 headername = new char[30] ;
289 strcpy(headername, "scalar_section.dxhead") ;
290 }
291 else {
292 headername = new char[ strlen(filename0)+9 ] ;
293 strcpy(headername, filename0) ;
294 strcat(headername, ".dxhead") ;
295 }
296
297 ofstream fheader(headername) ;
298
299 fheader << "file = " << filename << endl ;
300 fheader << "grid = " << nu << " x " << nv << endl ;
301 fheader << "format = ascii" << endl ;
302 fheader << "interleaving = record" << endl ;
303 fheader << "majority = column" << endl ;
304 fheader << "header = lines 4" << endl ;
305 fheader << "field = " << title_quotes << endl ;
306 fheader << "structure = scalar" << endl ;
307 fheader << "type = float" << endl ;
308 fheader << "dependency = positions" << endl ;
309 fheader << "positions = regular, regular, " << umin << ", " << du
310 << ", " << vmin << ", " << dv << endl ;
311 fheader << endl ;
312 fheader << "end" << endl ;
313
314 fheader.close() ;
315
316
317 if ( start_dx ) { // Launch of OpenDX
318
319 char* commande = new char[ strlen(headername) + 60 ] ;
320 strcpy(commande, "ln -s ") ;
321 strcat(commande, headername) ;
322 strcat(commande, " visu_section.dxhead") ;
323
324 system("rm -f visu_section.dxhead") ;
325 system(commande) ; // ln -s headername visu_section.general
326 system("dx -image visu_section.net &") ;
327
328 delete [] commande ;
329 }
330
331 // Final cleaning
332 // --------------
333 delete [] title ;
334 delete [] title_quotes ;
335 delete [] filename ;
336 delete [] headername ;
337
338}
339
340
341 //------------------------------//
342 // visu_box //
343 //------------------------------//
344
345void Scalar::visu_box(double xmin, double xmax, double ymin, double ymax,
346 double zmin, double zmax, const char* title0, const char* filename0,
347 bool start_dx, int nx, int ny, int nz) const {
348
349 const Scalar* scal ;
350 Scalar* scal_tmp = 0x0 ;
351
352 // Decrease of dzpuis if dzpuis != 0
353 if ( !check_dzpuis(0) ) {
354 scal_tmp = new Scalar(*this) ;
355 scal_tmp->dec_dzpuis(dzpuis) ;
356 scal = scal_tmp ;
357 }
358 else{
359 scal = this ;
360 }
361
362 char* title ;
363 char* title_quotes ;
364 if (title0 == 0x0) {
365 title = new char[2] ;
366 strcpy(title, " ") ;
367
368 title_quotes = new char[4] ;
369 strcpy(title_quotes, "\" \"") ;
370 }
371 else {
372 title = new char[ strlen(title0)+1 ] ;
373 strcpy(title, title0) ;
374
375 title_quotes = new char[ strlen(title0)+3 ] ;
376 strcpy(title_quotes, "\"") ;
377 strcat(title_quotes, title0) ;
378 strcat(title_quotes, "\"") ;
379 }
380
381 // --------------------------------------------------------
382 // Data file for OpenDX
383 // --------------------------------------------------------
384
385 char* filename ;
386 if (filename0 == 0x0) {
387 filename = new char[30] ;
388 strcpy(filename, "scalar_box.dxdata") ;
389 }
390 else {
391 filename = new char[ strlen(filename0)+8 ] ;
392 strcpy(filename, filename0) ;
393 strcat(filename, ".dxdata") ;
394 }
395
396 ofstream fdata(filename) ; // output file
397
398 fdata << title << "\n" ;
399 fdata << "size : " << nx << " x " << ny << " x " << nz << "\n" ;
400 fdata << "x_min = " << xmin << " x_max = " << xmax << "\n" ;
401 fdata << "y_min = " << ymin << " y_max = " << ymax << "\n" ;
402 fdata << "z_min = " << zmin << " z_max = " << zmax << "\n" ;
403
404 // The spectral coefficients are required
405 const Valeur& val = scal->va ;
406 val.coef() ;
407 const Mtbl_cf& cva = *(val.c_cf) ;
408
409 // What follows assumes that the mapping is radial:
410 assert( dynamic_cast<const Map_radial*>(mp) != 0x0 ) ;
411
412 fdata.precision(5) ;
413 fdata.setf(ios::scientific) ;
414
415 // Loop on the points of the drawing box
416 // ---------------------------------------
417 double dx = (xmax - xmin) / double(nx-1) ;
418 double dy = (ymax - ymin) / double(ny-1) ;
419 double dz = (zmax - zmin) / double(nz-1) ;
420
421 int npoint = 0 ; // number of data points per line in the file
422
423 for (int k=0; k<nz; k++) {
424
425 double zz = zmin + dz * k ;
426
427 for (int j=0; j<ny; j++) {
428
429 double yy = ymin + dy * j ;
430
431 for (int i=0; i<nx; i++) {
432
433 double xx = xmin + dx * i ;
434
435 // Values of (r,theta,phi) corresponding to (xa,ya,za) :
436 double rr, th, ph ; // polar coordinates of the mapping associated
437 // to *this
438
439 mp->convert_absolute(xx, yy, zz, rr, th, ph) ;
440
441 // Values of (l,xi,theta',phi') corresponding to (r,theta,phi):
442 double xi ; int l ;
443
444 mp->val_lx(rr, th, ph, l, xi) ; // radial mapping assumption
445
446 // Field value at this point:
447
448 double ff = cva.val_point(l, xi, th, ph) ;
449
450 fdata.width(14) ;
451 fdata << ff ;
452 npoint++ ;
453
454 if (npoint == 9) {
455 fdata << "\n" ;
456 npoint = 0 ;
457 }
458
459 }
460 }
461
462 }
463
464 if (npoint != 0) fdata << "\n" ;
465
466 fdata.close() ;
467
468 // --------------------------------------------------------
469 // Header file for OpenDX
470 // --------------------------------------------------------
471
472 char* headername ;
473 if (filename0 == 0x0) {
474 headername = new char[30] ;
475 strcpy(headername, "scalar_box.dxhead") ;
476 }
477 else {
478 headername = new char[ strlen(filename0)+9 ] ;
479 strcpy(headername, filename0) ;
480 strcat(headername, ".dxhead") ;
481 }
482
483 ofstream fheader(headername) ;
484
485 fheader << "file = " << filename << endl ;
486 fheader << "grid = " << nx << " x " << ny << " x " << nz << endl ;
487 fheader << "format = ascii" << endl ;
488 fheader << "interleaving = record" << endl ;
489 fheader << "majority = column" << endl ;
490 fheader << "header = lines 5" << endl ;
491 fheader << "field = " << title_quotes << endl ;
492 fheader << "structure = scalar" << endl ;
493 fheader << "type = float" << endl ;
494 fheader << "dependency = positions" << endl ;
495 fheader << "positions = regular, regular, regular, "
496 << xmin << ", " << dx << ", "
497 << ymin << ", " << dy << ", "
498 << zmin << ", " << dz << endl ;
499 fheader << endl ;
500 fheader << "end" << endl ;
501
502 fheader.close() ;
503
504
505 if ( start_dx ) { // Launch of OpenDX
506
507 char* commande = new char[ strlen(headername) + 60 ] ;
508 strcpy(commande, "ln -s ") ;
509 strcat(commande, headername) ;
510 strcat(commande, " visu_scalar_box.dxhead") ;
511
512 system("rm -f visu_scalar_box.dxhead") ;
513 system(commande) ; // ln -s headername visu_section.general
514 system("dx -image visu_scalar_box.net &") ;
515
516 delete [] commande ;
517 }
518
519
520 // Final cleaning
521 // --------------
522
523 if (scal_tmp != 0x0) delete scal_tmp ;
524 delete [] title ;
525 delete [] title_quotes ;
526 delete [] filename ;
527 delete [] headername ;
528
529
530}
531
532
533 //-------------------------------------//
534 // visu_section_anim //
535 //-------------------------------------//
536
537
538void Scalar::visu_section_anim(const char section_type, double aa, double umin,
539 double umax, double vmin, double vmax, int jtime, double ,
540 int jgraph, const char* title, const char* filename_root, bool start_dx,
541 int nu, int nv) const {
542
543 if ( jtime % jgraph != 0 ) return ;
544
545 // Preparation of the name of output file
546 // --------------------------------------
547 int k = jtime / jgraph ;
548
549 char* filename ;
550 if (filename_root == 0x0) {
551 filename = new char[40] ;
552 strcpy(filename, "anim") ;
553 }
554 else {
555 filename = new char[ strlen(filename_root)+10 ] ;
556 strcpy(filename, filename_root) ;
557 }
558
559 char nomk[5] ;
560 sprintf(nomk, "%04d", k) ;
561 strcat(filename, nomk) ;
562
563 // Call to visu_section to create the output file
564 // ----------------------------------------------
565
566 visu_section(section_type, aa, umin, umax, vmin, vmax, title, filename,
567 false, nu, nv) ;
568
569 // Shall we start OpenDX ?
570 // ---------------------
571
572 if ( start_dx ) { // Launch of OpenDX
573
574 system("dx -edit anime.net &") ;
575
576 }
577
578 // Final cleaning
579 // --------------
580
581 delete [] filename ;
582
583}
584}
Lorene::Map_radial
Base class for pure radial mappings.
Definition map.h:1551
Lorene::Mtbl_cf
Coefficients storage for the multi-domain spectral method.
Definition mtbl_cf.h:196
Lorene::Mtbl_cf::val_point
double val_point(int l, double x, double theta, double phi) const
Computes the value of the field represented by *this at an arbitrary point, by means of the spectral ...
Definition mtbl_cf_val_point.C:118
Lorene::Scalar::visu_box
void visu_box(double xmin, double xmax, double ymin, double ymax, double zmin, double zmax, const char *title0=0x0, const char *filename0=0x0, bool start_dx=true, int nx=40, int ny=40, int nz=40) const
3D visualization (volume rendering) via OpenDX.
Definition scalar_visu.C:345
Lorene::Scalar::check_dzpuis
bool check_dzpuis(int dzi) const
Returns false if the last domain is compactified and *this is not zero in this domain and dzpuis is n...
Definition scalar.C:879
Lorene::Scalar::Scalar
Scalar(const Map &mpi)
Constructor from mapping.
Definition scalar.C:210
Lorene::Scalar::visu_section_anim
void visu_section_anim(const char section_type, double aa, double umin, double umax, double vmin, double vmax, int jtime, double ttime, int jgraph=1, const char *title=0x0, const char *filename_root=0x0, bool start_dx=false, int nu=200, int nv=200) const
3D visualization via time evolving plane section (animation).
Definition scalar_visu.C:538
Lorene::Scalar::visu_section
void visu_section(const char section_type, double aa, double umin, double umax, double vmin, double vmax, const char *title=0x0, const char *filename=0x0, bool start_dx=true, int nu=200, int nv=200) const
3D visualization via a plane section.
Definition scalar_visu.C:81
Lorene::Scalar::va
Valeur va
The numerical value of the Scalar.
Definition scalar.h:411
Lorene::Scalar::dec_dzpuis
virtual void dec_dzpuis(int dec=1)
Decreases by dec units the value of dzpuis and changes accordingly the values of the Scalar in the co...
Definition scalar_r_manip.C:421
Lorene::Scalar::dzpuis
int dzpuis
Power of r by which the quantity represented by this must be divided in the compactified external d...
Definition scalar.h:409
Lorene::Tbl
Basic array class.
Definition tbl.h:161
Lorene::Tbl::set_etat_qcq
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition tbl.C:364
Lorene::Tbl::set
double & set(int i)
Read/write of a particular element (index i) (1D case).
Definition tbl.h:281
Lorene::Valeur
Values and coefficients of a (real-value) function.
Definition valeur.h:297
Lorene::Valeur::c_cf
Mtbl_cf * c_cf
Coefficients of the spectral expansion of the function.
Definition valeur.h:312
Lorene::Valeur::coef
void coef() const
Computes the coeffcients of *this.
Definition valeur_coef.C:151
Lorene::Tensor::mp
const Map *const mp
Mapping on which the numerical values at the grid points are defined.
Definition tensor.h:301
Lorene
Lorene prototypes.
Definition app_hor.h:67
Lorene::xa
Coord xa
Absolute x coordinate.
Definition map.h:742
Lorene::za
Coord za
Absolute z coordinate.
Definition map.h:744
Lorene::ya
Coord ya
Absolute y coordinate.
Definition map.h:743