ASCOT5
Loading...
Searching...
No Matches
hdf5_marker.c
Go to the documentation of this file.
1
8#include <stdio.h>
9#include <stdlib.h>
10#include <string.h>
11#include <hdf5.h>
12#include <hdf5_hl.h>
13#include "hdf5_helpers.h"
14#include "../particle.h"
15#include "../math.h"
16#include "../print.h"
17#include "../consts.h"
18#include "hdf5_marker.h"
19
20#define MRKPATH
22int hdf5_marker_read_particle(hid_t f, int* nmrk, input_particle** p,
23 char* qid);
24int hdf5_marker_read_guiding_center(hid_t f, int* nmrk, input_particle** p,
25 char* qid);
26int hdf5_marker_read_field_line(hid_t f, int* nmrk, input_particle** p,
27 char* qid);
28
41int hdf5_marker_read(hid_t f, int *n, input_particle** p, char* qid) {
42
43 char path[256];
44 int err = 1;
45
46 /* Read data the QID corresponds to */
47
48 hdf5_gen_path("/marker/prt_XXXXXXXXXX", qid, path);
49 if(hdf5_find_group(f, path) == 0) {
50 if( !hdf5_marker_read_particle(f, n, p, qid) ) {
51 print_out(VERBOSE_IO,"\nLoaded %d particles.\n", *n);
52 err = 0;
53 }
54 }
55
56 hdf5_generate_qid_path("/marker/gc_XXXXXXXXXX", qid, path);
57 if(hdf5_find_group(f, path) == 0) {
58 if( !hdf5_marker_read_guiding_center(f, n, p, qid) ) {
59 print_out(VERBOSE_IO,"\nLoaded %d guiding centers.\n", *n);
60 err = 0;
61 }
62 }
63
64 hdf5_generate_qid_path("/marker/fl_XXXXXXXXXX", qid, path);
65 if(hdf5_find_group(f, path) == 0) {
66 if( !hdf5_marker_read_field_line(f, n, p, qid) ) {
67 print_out(VERBOSE_IO,"\nLoaded %d field lines.\n", *n);
68 err = 0;
69 }
70 }
71
72 return err;
73}
74
87int hdf5_marker_read_particle(hid_t f, int* nmrk, input_particle** mrk,
88 char* qid) {
89 #undef MRKPATH
90 #define MRKPATH "/marker/prt_XXXXXXXXXX/"
91
92 integer n;
93 if( hdf5_read_long(MRKPATH "n", &n,
94 f, qid, __FILE__, __LINE__) ) {return 1;}
95 *nmrk = n;
96
97 real* r = malloc(n * sizeof(real));
98 real* phi = malloc(n * sizeof(real));
99 real* z = malloc(n * sizeof(real));
100 real* v_r = malloc(n * sizeof(real));
101 real* v_phi = malloc(n * sizeof(real));
102 real* v_z = malloc(n * sizeof(real));
103 real* mass = malloc(n * sizeof(real));
104 int* charge = malloc(n * sizeof(int));
105 int* anum = malloc(n * sizeof(int));
106 int* znum = malloc(n * sizeof(int));
107 real* weight = malloc(n * sizeof(real));
108 real* time = malloc(n * sizeof(real));
109 integer* id = malloc(n * sizeof(integer));
110
111 if( hdf5_read_double(MRKPATH "r", r,
112 f, qid, __FILE__, __LINE__) ) {return 1;}
113 if( hdf5_read_double(MRKPATH "phi", phi,
114 f, qid, __FILE__, __LINE__) ) {return 1;}
115 if( hdf5_read_double(MRKPATH "z", z,
116 f, qid, __FILE__, __LINE__) ) {return 1;}
117 if( hdf5_read_double(MRKPATH "vr", v_r,
118 f, qid, __FILE__, __LINE__) ) {return 1;}
119 if( hdf5_read_double(MRKPATH "vphi", v_phi,
120 f, qid, __FILE__, __LINE__) ) {return 1;}
121 if( hdf5_read_double(MRKPATH "vz", v_z,
122 f, qid, __FILE__, __LINE__) ) {return 1;}
123 if( hdf5_read_double(MRKPATH "mass", mass,
124 f, qid, __FILE__, __LINE__) ) {return 1;}
125 if( hdf5_read_int(MRKPATH "charge", charge,
126 f, qid, __FILE__, __LINE__) ) {return 1;}
127 if( hdf5_read_int(MRKPATH "anum", anum,
128 f, qid, __FILE__, __LINE__) ) {return 1;}
129 if( hdf5_read_int(MRKPATH "znum", znum,
130 f, qid, __FILE__, __LINE__) ) {return 1;}
131 if( hdf5_read_double(MRKPATH "weight", weight,
132 f, qid, __FILE__, __LINE__) ) {return 1;}
133 if( hdf5_read_double(MRKPATH "time", time,
134 f, qid, __FILE__, __LINE__) ) {return 1;}
135 if( hdf5_read_long(MRKPATH "id", id,
136 f, qid, __FILE__, __LINE__) ) {return 1;}
137
138 *mrk = (input_particle*) malloc(n * sizeof(input_particle));
139 input_particle* p = *mrk;
140
141 for(integer i = 0; i < n; i++) {
142 real gamma = 1.0/sqrt(
143 1 - (v_r[i]*v_r[i] + v_phi[i]*v_phi[i] + v_z[i]*v_z[i]) / CONST_C2 );
144 p[i].p.r = r[i];
145 p[i].p.phi = phi[i] * CONST_PI / 180;
146 p[i].p.z = z[i];
147 p[i].p.p_r = v_r[i] * gamma * mass[i]* CONST_U;;
148 p[i].p.p_phi = v_phi[i] * gamma * mass[i]* CONST_U;;
149 p[i].p.p_z = v_z[i] * gamma * mass[i]* CONST_U;;
150 p[i].p.mass = mass[i] * CONST_U;
151 p[i].p.charge = charge[i] * CONST_E;
152 p[i].p.anum = anum[i];
153 p[i].p.znum = znum[i];
154 p[i].p.weight = weight[i];
155 p[i].p.time = time[i];
156 p[i].p.id = (integer) id[i];
158 }
159
160 free(r);
161 free(phi);
162 free(z);
163 free(v_r);
164 free(v_phi);
165 free(v_z);
166 free(mass);
167 free(charge);
168 free(anum);
169 free(znum);
170 free(weight);
171 free(time);
172 free(id);
173
174 return 0;
175}
176
190 char* qid) {
191 #undef MRKPATH
192 #define MRKPATH "/marker/gc_XXXXXXXXXX/"
193
194 integer n;
195 if( hdf5_read_long(MRKPATH "n", &n,
196 f, qid, __FILE__, __LINE__) ) {return 1;}
197 *nmrk = n;
198
199 real* r = malloc(n * sizeof(real));
200 real* phi = malloc(n * sizeof(real));
201 real* z = malloc(n * sizeof(real));
202 real* energy = malloc(n * sizeof(real));
203 real* pitch = malloc(n * sizeof(real));
204 real* zeta = malloc(n * sizeof(real));
205 real* mass = malloc(n * sizeof(real));
206 int* charge = malloc(n * sizeof(int));
207 int* anum = malloc(n * sizeof(int));
208 int* znum = malloc(n * sizeof(int));
209 real* weight = malloc(n * sizeof(real));
210 real* time = malloc(n * sizeof(real));
211 integer* id = malloc(n * sizeof(integer));
212
213 if( hdf5_read_double(MRKPATH "r", r,
214 f, qid, __FILE__, __LINE__) ) {return 1;}
215 if( hdf5_read_double(MRKPATH "phi", phi,
216 f, qid, __FILE__, __LINE__) ) {return 1;}
217 if( hdf5_read_double(MRKPATH "z", z,
218 f, qid, __FILE__, __LINE__) ) {return 1;}
219 if( hdf5_read_double(MRKPATH "energy", energy,
220 f, qid, __FILE__, __LINE__) ) {return 1;}
221 if( hdf5_read_double(MRKPATH "pitch", pitch,
222 f, qid, __FILE__, __LINE__) ) {return 1;}
223 if( hdf5_read_double(MRKPATH "zeta", zeta,
224 f, qid, __FILE__, __LINE__) ) {return 1;}
225 if( hdf5_read_double(MRKPATH "mass", mass,
226 f, qid, __FILE__, __LINE__) ) {return 1;}
227 if( hdf5_read_int(MRKPATH "charge", charge,
228 f, qid, __FILE__, __LINE__) ) {return 1;}
229 if( hdf5_read_int(MRKPATH "anum", anum,
230 f, qid, __FILE__, __LINE__) ) {return 1;}
231 if( hdf5_read_int(MRKPATH "znum", znum,
232 f, qid, __FILE__, __LINE__) ) {return 1;}
233 if( hdf5_read_double(MRKPATH "weight", weight,
234 f, qid, __FILE__, __LINE__) ) {return 1;}
235 if( hdf5_read_double(MRKPATH "time", time,
236 f, qid, __FILE__, __LINE__) ) {return 1;}
237 if( hdf5_read_long(MRKPATH "id", id,
238 f, qid, __FILE__, __LINE__) ) {return 1;}
239
240 *mrk = (input_particle*) malloc(n * sizeof(input_particle));
241 input_particle* p = *mrk;
242
243 for(integer i = 0; i < n; i++) {
244 p[i].p_gc.r = r[i];
245 p[i].p_gc.phi = phi[i] * CONST_PI / 180;
246 p[i].p_gc.z = z[i];
247 p[i].p_gc.energy = energy[i] * CONST_E;
248 p[i].p_gc.pitch = pitch[i];
249 p[i].p_gc.zeta = zeta[i];
250 p[i].p_gc.mass = mass[i] * CONST_U;
251 p[i].p_gc.charge = charge[i] * CONST_E;
252 p[i].p_gc.anum = anum[i];
253 p[i].p_gc.znum = znum[i];
254 p[i].p_gc.weight = weight[i];
255 p[i].p_gc.time = time[i];
256 p[i].p_gc.id = (integer) id[i];
258 }
259
260 free(r);
261 free(phi);
262 free(z);
263 free(energy);
264 free(pitch);
265 free(zeta);
266 free(mass);
267 free(charge);
268 free(anum);
269 free(znum);
270 free(weight);
271 free(time);
272 free(id);
273
274 return 0;
275}
276
289int hdf5_marker_read_field_line(hid_t f, int* nmrk, input_particle** mrk,
290 char* qid) {
291 #undef MRKPATH
292 #define MRKPATH "/marker/fl_XXXXXXXXXX/"
293
294 integer n;
295 if( hdf5_read_long(MRKPATH "n", &n,
296 f, qid, __FILE__, __LINE__) ) {return 1;}
297 *nmrk = n;
298
299 real* r = malloc(n * sizeof(real));
300 real* phi = malloc(n * sizeof(real));
301 real* z = malloc(n * sizeof(real));
302 real* pitch = malloc(n * sizeof(real));
303 real* weight = malloc(n * sizeof(real));
304 real* time = malloc(n * sizeof(real));
305 integer* id = malloc(n * sizeof(integer));
306
307 if( hdf5_read_double(MRKPATH "r", r,
308 f, qid, __FILE__, __LINE__) ) {return 1;}
309 if( hdf5_read_double(MRKPATH "phi", phi,
310 f, qid, __FILE__, __LINE__) ) {return 1;}
311 if( hdf5_read_double(MRKPATH "z", z,
312 f, qid, __FILE__, __LINE__) ) {return 1;}
313 if( hdf5_read_double(MRKPATH "pitch", pitch,
314 f, qid, __FILE__, __LINE__) ) {return 1;}
315 if( hdf5_read_double(MRKPATH "weight", weight,
316 f, qid, __FILE__, __LINE__) ) {return 1;}
317 if( hdf5_read_double(MRKPATH "time", time,
318 f, qid, __FILE__, __LINE__) ) {return 1;}
319 if( hdf5_read_long(MRKPATH "id", id,
320 f, qid, __FILE__, __LINE__) ) {return 1;}
321
322 *mrk = (input_particle*) malloc(n * sizeof(input_particle));
323 input_particle* p = *mrk;
324
325 for(integer i = 0; i < n; i++) {
326 p[i].p_ml.r = r[i];
327 p[i].p_ml.phi = phi[i] * CONST_PI / 180;
328 p[i].p_ml.z = z[i];
329 p[i].p_ml.pitch = pitch[i];
330 p[i].p_ml.weight = weight[i];
331 p[i].p_ml.time = time[i];
332 p[i].p_ml.id = (integer)id[i];
334 }
335
336 free(r);
337 free(phi);
338 free(z);
339 free(pitch);
340 free(weight);
341 free(time);
342 free(id);
343
344 return 0;
345}
346
359int hdf5_marker_write_particle(hid_t f, int n, input_particle* p, char* qid) {
360
361 if(hdf5_find_group(f, "/marker/")) {
362 hdf5_create_group(f, "/marker/");
363 }
364
365 char path[256];
366 hdf5_gen_path("/marker/prt_XXXXXXXXXX", qid, path);
367
368 hid_t grp = H5Gcreate2(f, path, H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
369
370 real* r = malloc(n * sizeof(real));
371 real* phi = malloc(n * sizeof(real));
372 real* z = malloc(n * sizeof(real));
373 real* v_r = malloc(n * sizeof(real));
374 real* v_phi = malloc(n * sizeof(real));
375 real* v_z = malloc(n * sizeof(real));
376 real* mass = malloc(n * sizeof(real));
377 int* charge = malloc(n * sizeof(real));
378 int* anum = malloc(n * sizeof(int));
379 int* znum = malloc(n * sizeof(int));
380 real* weight = malloc(n * sizeof(real));
381 real* time = malloc(n * sizeof(real));
382 integer* id = malloc(n * sizeof(integer));
383
384 for(int i = 0; i < n; i++) {
385 real gamma = sqrt(
386 1
387 + (p[i].p.p_r*p[i].p.p_r
388 + p[i].p.p_phi*p[i].p.p_phi
389 + p[i].p.p_z*p[i].p.p_z) / (CONST_C2*p[i].p.mass)
390 );
391 r[i] = p[i].p.r;
392 phi[i] = p[i].p.phi * 180 / CONST_PI;
393 z[i] = p[i].p.z;
394 v_r[i] = p[i].p.p_r / (p[i].p.mass*gamma);
395 v_phi[i] = p[i].p.p_phi / (p[i].p.mass*gamma);
396 v_z[i] = p[i].p.p_z / (p[i].p.mass*gamma);
397 mass[i] = p[i].p.mass / CONST_U;
398 charge[i] = round(p[i].p.charge / CONST_E);
399 anum[i] = p[i].p.anum;
400 znum[i] = p[i].p.znum;
401 weight[i] = p[i].p.weight;
402 time[i] = p[i].p.time;
403 id[i] = p[i].p.id;
404 }
405
406 hdf5_write_extendible_dataset_int(grp, "n", 1, &n);
408 hdf5_write_extendible_dataset_double(grp, "phi", n, phi);
410 hdf5_write_extendible_dataset_double(grp, "vr", n, v_r);
411 hdf5_write_extendible_dataset_double(grp, "vphi", n, v_phi);
412 hdf5_write_extendible_dataset_double(grp, "vz", n, v_z);
413 hdf5_write_extendible_dataset_double(grp, "mass", n, mass);
414 hdf5_write_extendible_dataset_int(grp, "charge", n, charge);
415 hdf5_write_extendible_dataset_int(grp, "anum", n, anum);
416 hdf5_write_extendible_dataset_int(grp, "znum", n, znum);
417 hdf5_write_extendible_dataset_double(grp, "weight", n, weight);
418 hdf5_write_extendible_dataset_double(grp, "time", n, time);
419 hdf5_write_extendible_dataset_long(grp, "id", n, id);
420
421 H5Gclose(grp);
422
423 free(r);
424 free(phi);
425 free(z);
426 free(v_r);
427 free(v_phi);
428 free(v_z);
429 free(mass);
430 free(charge);
431 free(anum);
432 free(znum);
433 free(weight);
434 free(time);
435 free(id);
436
437 return 0;
438}
double real
Definition ascot5.h:85
long integer
Definition ascot5.h:84
Header file containing physical and mathematical constants.
#define CONST_U
Atomic mass unit in kilograms [kg]
Definition consts.h:29
#define CONST_PI
pi
Definition consts.h:11
#define CONST_E
Elementary charge [C]
Definition consts.h:32
#define CONST_C2
Speed of light squared [m^2/s^2]
Definition consts.h:26
hid_t hdf5_create_group(hid_t loc, const char *path)
Create a group (with parent groups if necessary). Returns a handle to the group. Negative on failure.
herr_t hdf5_find_group(hid_t loc, const char *path)
Checks if given group exists within given hdf5 file. Negative value is returned if the group doesn't ...
herr_t hdf5_write_extendible_dataset_long(hid_t group, const char *datasetname, int length, long *data)
Create and write to an extendible dataset for long data.
herr_t hdf5_write_extendible_dataset_int(hid_t group, const char *datasetname, int length, int *data)
Create and write to an extendible dataset int data.
int hdf5_read_double(const char *var, real *ptr, hid_t file, char *qid, const char *errfile, int errline)
Read double-valued data from ASCOT5 HDF5 file.
int hdf5_read_long(const char *var, long *ptr, hid_t file, char *qid, const char *errfile, int errline)
Read long-valued data from ASCOT5 HDF5 file.
char * hdf5_gen_path(const char *original, char *qid, char *path)
Generate a valid path from a given template and qid.
herr_t hdf5_write_extendible_dataset_double(hid_t group, const char *datasetname, int length, double *data)
Create and write to an extendible dataset for double data.
int hdf5_read_int(const char *var, int *ptr, hid_t file, char *qid, const char *errfile, int errline)
Read int-valued data from ASCOT5 HDF5 file.
char * hdf5_generate_qid_path(const char *original, char *qid, char *path)
Generate a valid path from a given template and qid.
Header file for hdf5_helpers.h.
int hdf5_marker_read_field_line(hid_t f, int *nmrk, input_particle **p, char *qid)
Read field line input.
int hdf5_marker_read_guiding_center(hid_t f, int *nmrk, input_particle **p, char *qid)
Read guiding center input.
int hdf5_marker_read(hid_t f, int *n, input_particle **p, char *qid)
Read marker input.
Definition hdf5_marker.c:41
#define MRKPATH
Definition hdf5_marker.c:20
int hdf5_marker_write_particle(hid_t f, int n, input_particle *p, char *qid)
Write particle input.
int hdf5_marker_read_particle(hid_t f, int *nmrk, input_particle **p, char *qid)
Read particle input.
Definition hdf5_marker.c:87
Header file for hdf5_marker.c.
Header file for math.c.
Header file for particle.c.
@ input_particle_type_gc
Definition particle.h:170
@ input_particle_type_p
Definition particle.h:169
@ input_particle_type_ml
Definition particle.h:171
Macros for printing console output.
#define print_out(v,...)
Print to standard output.
Definition print.h:31
@ VERBOSE_IO
Definition print.h:20
Wrapper for marker structs.
Definition particle.h:186
particle_ml p_ml
Definition particle.h:191
input_particle_type type
Definition particle.h:187
particle p
Definition particle.h:189
particle_gc p_gc
Definition particle.h:190
real weight
Definition particle.h:121
real energy
Definition particle.h:114
integer id
Definition particle.h:123
real charge
Definition particle.h:118
real weight
Definition particle.h:137
integer id
Definition particle.h:139
real charge
Definition particle.h:96
int anum
Definition particle.h:97
real z
Definition particle.h:91
real weight
Definition particle.h:99
real phi
Definition particle.h:90
real p_phi
Definition particle.h:93
real time
Definition particle.h:100
int znum
Definition particle.h:98
real r
Definition particle.h:89
real mass
Definition particle.h:95
real p_z
Definition particle.h:94
integer id
Definition particle.h:101
real p_r
Definition particle.h:92