modwritemeta.f90 20.4 KB
Newer Older
1
2
module modwritemeta

3
  use iso_fortran_env, only : OUTPUT_UNIT
4

5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
  use modconstant, only : type_parameter_halofinder,&
       type_parameter_pfof,&
       type_parameter_pfof_snap,&
       type_parameter_pfof_cone,&
       type_parameter_conecreator,&
       type_parameter_conecreator_part,&
       type_parameter_conecreator_grav,&
       type_parameter_psod_snap,&
       type_info_ramses,&
       type_info_cone,&
       type_info_cone_part,&
       type_info_cone_grav,&
       type_common_metadata,&
       NAME_PFOF_SNAP,&
       NAME_PFOF_CONE,&
       NAME_PSOD_SNAP,&
       NAME_CONECREATOR_PART,&
       NAME_CONECREATOR_GRAV,&
23
       GIT_VERSION
24

25
  use modhdf5, only : hid_t,&
26
27
28
29
30
31
32
33
       H5_STR_LEN,&
       Hdf5_create_group,&
       Hdf5_open_group,&
       Hdf5_close_group,&
       Hdf5_read_attr,&
       Hdf5_read_data,&
       Hdf5_write_attr,&
       Hdf5_write_data
34

35
36
  use modmpicommons, only : procid

37
38
39
  implicit none

  private
40
41
42
43
44
  public :: Write_meta_common, &
       Write_meta_halofinder_parameter, &
       Write_meta_conecreator_parameter, &
       Write_meta_info_ramses, &
       Write_meta_info_cone
45
46

  ! Some metadata constants whose value may change for some software or options
47
48
  character(len=H5_STR_LEN) :: PARTICLE_TYPE = 'dark_matter'
  character(len=H5_STR_LEN) :: PFOF_CELL_ORDER = 'none'
49
50
51
52

contains
  !=======================================================================
  !> Write metadata common to every hdf5 ouput files of the pfof toolbox
53
  subroutine Write_meta_common(file_id, codename, npart, constant_mass, process_id)
54
55

    integer(kind=hid_t), intent(in) :: file_id
56
    character(len=H5_STR_LEN), intent(in) :: codename
57
    integer(kind=8), intent(in) :: npart
58
    integer(kind=4), intent(in), optional :: constant_mass
59
    integer(kind=4), intent(in), optional :: process_id
60

61
62
63
64
    character(len=H5_STR_LEN) :: groupname
    character(len=H5_STR_LEN) :: aname
    character(len=H5_STR_LEN) :: adata
    character(len=H5_STR_LEN) :: dsetname
65
66
67
    integer(kind=hid_t) :: gr_id
    integer(kind=4) :: tmpint4

68
69
70
#ifdef DEBUG
    write(OUTPUT_UNIT,*) 'Enter Write_meta_common on process ', procid
#endif
71

72
    groupname = 'metadata'
73
    call Hdf5_open_group(file_id,groupname, gr_id)
74
75
    ! Common metadata
    aname = 'created_by'
76
    call Hdf5_write_attr(gr_id, aname, codename)
77
78
    aname = 'version'
    call Hdf5_write_attr(gr_id, aname, GIT_VERSION)
79
80
    aname = 'simulation_code'
    adata = 'ramses'
81
    call Hdf5_write_attr(gr_id, aname,adata)
82
83
    aname = 'particle_type'
    adata = PARTICLE_TYPE
84
    call Hdf5_write_attr(gr_id, aname, adata)
85
86
87
88
89
90
    aname = 'process_id'
    if(present(process_id)) then
       tmpint4=process_id
    else
       tmpint4=0
    end if
91
    call Hdf5_write_attr(gr_id, aname, tmpint4)
92
    aname = 'constant_mass'
93
94
95
96
97
    if(present(constant_mass)) then
       tmpint4 = constant_mass
    else
       tmpint4 = -1
    end if
98
    call Hdf5_write_attr(gr_id, aname, tmpint4)
99
100
    aname = 'units'
    adata = 'ramses'
101
    call Hdf5_write_attr(gr_id, aname, adata)
102
103
    aname = 'pfof_cells_order'
    adata = PFOF_CELL_ORDER
104
    call Hdf5_write_attr(gr_id, aname, adata)
105
106
    ! Write the number of particles in this simulation as an integer(kind=8) dataset
    dsetname = 'npart_simulation'
107
108
    call Hdf5_write_data(gr_id, dsetname, npart)
    call Hdf5_close_group(gr_id)
109

110
111
112
113
114
#ifdef DEBUG
    write(OUTPUT_UNIT,*) 'Exit Write_meta_common on process ', procid
#endif

  end subroutine Write_meta_common
115
116

  !=======================================================================
117
  !> Write halofinder input parameters as metadata
118
  subroutine Write_meta_halofinder_parameter(file_id, param_halofinder)
119
120
121
122

    integer(kind=hid_t), intent(in) :: file_id
    class(type_parameter_halofinder), intent(in) :: param_halofinder

123
124
    character(len=H5_STR_LEN) :: groupname
    character(len=H5_STR_LEN) :: aname
125
126
127
128
129
130
131
132
    integer(kind=hid_t) :: gr_id
    integer(kind=hid_t) :: gr_pfof_id
    integer(kind=hid_t) :: gr_input_id
    integer(kind=hid_t) :: gr_fof_id
    integer(kind=hid_t) :: gr_output_id
    integer(kind=4) :: tmpint4
    character(len=16) :: c_b_name

133
134
135
136
#ifdef DEBUG
    write(OUTPUT_UNIT,*) 'Enter Write_meta_halofinder_parameter on process ', procid
#endif

137
138
139
    select type (param_halofinder)
    type is (type_parameter_pfof_snap)
       c_b_name = NAME_PFOF_SNAP
140
    type is (type_parameter_pfof_cone)
141
       c_b_name = NAME_PFOF_CONE
142
    type is (type_parameter_psod_snap)
143
       c_b_name = NAME_PSOD_SNAP
144
    end select
145
    groupname = 'metadata'
146
    call Hdf5_open_group(file_id,groupname, gr_id)
147
148
    ! pfof parameters:
    groupname = trim(c_b_name)//'_parameters'
149
    call Hdf5_create_group(gr_id, groupname, gr_pfof_id)
150
    groupname = 'input_parameters'
151
    call Hdf5_create_group(gr_pfof_id, groupname, gr_input_id)
152
    aname = 'input_path'
153
    call Hdf5_write_attr(gr_input_id, aname, param_halofinder%input_path)
154
    aname = 'part_input_file'
155
    call Hdf5_write_attr(gr_input_id, aname, param_halofinder%part_input_file)
156
157
158
159
160
161
162
    ! Write potential logical as an integer attribute (1=true, 0=false)
    aname = 'do_read_potential'
    if(param_halofinder%do_read_potential) then
       tmpint4=1
    else
       tmpint4=0
    end if
163
    call Hdf5_write_attr(gr_input_id, aname, tmpint4)
164
165
    ! Write force logical as an integer attribute (1=true, 0=false)
    aname = 'do_read_gravitational_field'
166
    if(param_halofinder%do_read_gravitational_field) then
167
       tmpint4=1
168
    else
169
       tmpint4=0
170
171
172
173
    end if
    call Hdf5_write_attr(gr_input_id, aname, tmpint4)
    select type (param_halofinder)
    type is (type_parameter_pfof_snap)
174
       aname = 'info_input_file'
175
       call Hdf5_write_attr(gr_input_id, aname, param_halofinder%info_input_file)
176
       aname = 'grpsize'
177
       call Hdf5_write_attr(gr_input_id, aname, param_halofinder%grpsize)
178
       aname = 'code_index'
179
       call Hdf5_write_attr(gr_input_id, aname, param_halofinder%code_index)
180
       aname = 'do_skip_star'
181
       if(param_halofinder%do_skip_star) then
182
          tmpint4=1
183
       else
184
          tmpint4=0
185
186
       end if
       call Hdf5_write_attr(gr_input_id, aname, tmpint4)
187
       aname = 'do_skip_metal'
188
       if(param_halofinder%do_skip_metal) then
189
          tmpint4=1
190
       else
191
          tmpint4=0
192
       end if
193
       call Hdf5_write_attr(gr_input_id, aname, tmpint4)
194
195
       aname='do_read_from_cube'
       tmpint4 = 0
196
197
       if(param_halofinder%do_read_from_cube) tmpint4 = 1
       call Hdf5_write_attr(gr_input_id, aname, tmpint4)
198
       aname = 'gatherread_factor'
199
200
       call Hdf5_write_attr(gr_input_id, aname, param_halofinder%gatherread_factor)
    type is (type_parameter_pfof_cone)
201
       aname = 'shell_first_id'
202
       call Hdf5_write_attr(gr_input_id, aname, param_halofinder%shell_first_id)
203
       aname = 'shell_last_id'
204
205
206
       call Hdf5_write_attr(gr_input_id, aname, param_halofinder%shell_last_id)
    end select
    call Hdf5_close_group(gr_input_id)
207
    groupname = 'fof_parameters'
208
    call Hdf5_create_group(gr_pfof_id, groupname, gr_fof_id)
209
210
    ! Write minimum halo mass Mmin as attribute
    aname = 'npart_halo_min'
211
    call Hdf5_write_attr(gr_fof_id, aname, param_halofinder%mmin)
212
    aname = 'npart_halo_max'
213
    call Hdf5_write_attr(gr_fof_id, aname, param_halofinder%mmax)
214
    ! Write doUnbinding as attribute (not implemented yet)
215
    aname = 'do_unbinding'
216
    if(param_halofinder%do_unbinding) then
217
       tmpint4=1
218
    else
219
       tmpint4=0
220
221
    end if
    call Hdf5_write_attr(gr_fof_id, aname, tmpint4)
222
    ! Write doSubHalo as attribute (not implemented yet)
223
    aname = 'do_subhalo'
224
    if(param_halofinder%do_subhalo) then
225
       tmpint4=1
226
    else
227
       tmpint4=0
228
229
230
231
    end if
    call Hdf5_write_attr(gr_fof_id, aname, tmpint4)
    select type (param_halofinder)
    type is (type_parameter_pfof_snap)
232
233
       ! Write percolation parameter b as attribute
       aname = 'percolation_length'
234
       call Hdf5_write_attr(gr_fof_id, aname, param_halofinder%percolation_length)
235
       aname = 'do_fof'
236
       if(param_halofinder%do_fof) then
237
          tmpint4=1
238
       else
239
          tmpint4=0
240
241
242
       end if
       call Hdf5_write_attr(gr_fof_id, aname, tmpint4)
    type is (type_parameter_pfof_cone)
243
244
       ! Write percolation parameter b as attribute
       aname = 'percolation_length'
245
       call Hdf5_write_attr(gr_fof_id, aname, param_halofinder%percolation_length)
246
    type is (type_parameter_psod_snap)
247
       aname = 'do_sod'
248
       if(param_halofinder%do_sod) then
249
          tmpint4=1
250
       else
251
          tmpint4=0
252
253
254
       end if
       call Hdf5_write_attr(gr_fof_id, aname, tmpint4)
    end select
255
    call Hdf5_close_group(gr_fof_id)
256
    groupname = 'output_parameters'
257
    call Hdf5_create_group(gr_pfof_id, groupname, gr_output_id)
258
259
    ! Write the simulation name as attribute
    aname = 'simulation_name'
260
    call Hdf5_write_attr(gr_output_id, aname, param_halofinder%simulation_name)
261
262
    ! Write do_timings as attribute
    aname = 'do_timings'
263
    if(param_halofinder%do_timings) then
264
       tmpint4=1
265
    else
266
       tmpint4=0
267
268
269
270
    end if
    call Hdf5_write_attr(gr_output_id, aname, tmpint4)
    select type (param_halofinder)
    type is (type_parameter_pfof_snap)
271
272
       ! write the snapshot number: if halo computed from a lightcone then outputNB = -1
       aname = 'snapshot'
273
       call Hdf5_write_attr(gr_output_id, aname, param_halofinder%snapshot)
274
       tmpint4 = 0
275
       if(param_halofinder%do_write_cube) tmpint4 = 1
276
       aname = 'do_write_cube'
277
       call Hdf5_write_attr(gr_output_id, aname, tmpint4)
278
       tmpint4 = 0
279
       if(param_halofinder%do_sort_cube) tmpint4 = 1
280
       aname = 'do_sort_cube'
281
       call Hdf5_write_attr(gr_output_id, aname, tmpint4)
282
       aname = 'gatherwrite_factor'
283
284
285
286
287
       call Hdf5_write_attr(gr_output_id, aname, param_halofinder%gatherwrite_factor)
    end select
    call Hdf5_close_group(gr_output_id)
    call Hdf5_close_group(gr_pfof_id)
    call Hdf5_close_group(gr_id)
288

289
290
291
#ifdef DEBUG
    write(OUTPUT_UNIT,*) 'Exit Write_meta_halofinder_parameter on process ', procid
#endif
292

293
  end subroutine Write_meta_halofinder_parameter
294
295
296
297


  !=======================================================================
  !> Write conecreator input parameters as metadata
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
  subroutine Write_meta_conecreator_parameter(file_id, param_cone)

    integer(kind=hid_t), intent(in) :: file_id
    class(type_parameter_conecreator), intent(in) :: param_cone

    character(len=H5_STR_LEN) :: groupname
    character(len=H5_STR_LEN) :: aname
    integer(kind=hid_t) :: gr_id
    integer(kind=hid_t) :: gr_param_id
    integer(kind=hid_t) :: gr_input_id
    integer(kind=hid_t) :: gr_output_id
    integer(kind=4) :: tmpint4
    character(len=16) :: c_b_name

#ifdef DEBUG
    write(OUTPUT_UNIT,*) 'Enter Write_meta_conecreator_parameter on process ', procid
#endif
315

316
317
    select type(param_cone)
    type is(type_parameter_conecreator_part)
318
       c_b_name = NAME_CONECREATOR_PART
319
    type is(type_parameter_conecreator_grav)
320
       c_b_name = NAME_CONECREATOR_GRAV
321
    end select
322
    groupname = 'metadata'
323
    call Hdf5_open_group(file_id, groupname, gr_id)
324
    groupname = c_b_name//'_parameters'
325
    call Hdf5_create_group(gr_id, groupname, gr_param_id)
326
    groupname = 'input_parameters'
327
    call Hdf5_create_group(gr_param_id, groupname, gr_input_id)
328
    aname = 'input_path'
329
    call Hdf5_write_attr(gr_input_id, aname, param_cone%input_path)
330
    aname = 'cone_input_file'
331
    call Hdf5_write_attr(gr_input_id, aname, param_cone%cone_input_file)
332
    aname = 'info_cone_input_file'
333
    call Hdf5_write_attr(gr_input_id, aname, param_cone%info_cone_input_file)
334
    aname = 'info_ramses_input_file'
335
    call Hdf5_write_attr(gr_input_id, aname, param_cone%info_ramses_input_file)
336
    aname = 'nfile'
337
    call Hdf5_write_attr(gr_input_id, aname, param_cone%nfile)
338
    aname = 'first_file'
339
    call Hdf5_write_attr(gr_input_id, aname, param_cone%first_file)
340
    aname = 'last_file'
341
    call Hdf5_write_attr(gr_input_id, aname, param_cone%last_file)
342
    aname = 'cone_max_radius'
343
344
345
    call Hdf5_write_attr(gr_input_id, aname, param_cone%cone_max_radius)
    select type(param_cone)
    type is (type_parameter_conecreator_part)
346
       aname = 'do_read_ramses_part_id'
347
       tmpint4 = 0
348
349
       if(param_cone%do_read_ramses_part_id) tmpint4 = 1
       call Hdf5_write_attr(gr_input_id, aname, tmpint4)
350
       aname = 'do_read_potential'
351
       tmpint4 = 0
352
353
       if(param_cone%do_read_potential) tmpint4 = 1
       call Hdf5_write_attr(gr_input_id, aname, tmpint4)
354
       aname = 'do_read_gravitational_field'
355
       tmpint4 = 0
356
357
358
       if(param_cone%do_read_gravitational_field) tmpint4 = 1
       call Hdf5_write_attr(gr_input_id, aname, tmpint4)
    type is (type_parameter_conecreator_grav)
359
       aname = 'nlevel'
360
       call Hdf5_write_attr(gr_input_id, aname, param_cone%nlevel)
361
       aname = 'levelmin'
362
363
364
       call Hdf5_write_attr(gr_input_id, aname, param_cone%levelmin)
    end select
    call Hdf5_close_group(gr_input_id)
365
    groupname = 'output_parameters'
366
    call Hdf5_create_group(gr_param_id, groupname, gr_output_id)
367
    aname = 'simulation_name'
368
    call Hdf5_write_attr(gr_output_id, aname, param_cone%simulation_name)
369
    aname = 'cube_size'
370
371
372
373
    call Hdf5_write_attr(gr_output_id, aname, param_cone%cube_size)
    call Hdf5_close_group(gr_output_id)
    call Hdf5_close_group(gr_param_id)
    call Hdf5_close_group(gr_id)
374

375
376
377
#ifdef DEBUG
    write(OUTPUT_UNIT,*) 'Exit Write_meta_conecreator_parameter on process ', procid
#endif
378

379
  end subroutine Write_meta_conecreator_parameter
380
381
382

  !=======================================================================
  !> Write ramses info as metadata
383
384
385
386
387
  subroutine Write_meta_info_ramses(file_id, inforamses, islast)

    integer(kind=hid_t), intent(in) :: file_id
    type(type_info_ramses), intent(in) :: inforamses
    logical(kind=4), intent(in) :: islast
388

389
390
391
392
    character(len=H5_STR_LEN) :: groupname
    character(len=H5_STR_LEN) :: aname
    integer(kind=hid_t) :: gr_id
    integer(kind=hid_t) :: gr_ramses_id
393

394
395
396
#ifdef DEBUG
    write(OUTPUT_UNIT,*) 'Enter Write_meta_info_ramses on process ', procid
#endif
397
398

    groupname = 'metadata'
399
    call Hdf5_open_group(file_id, groupname, gr_id)
400
401
    ! Ramses Info Metadata
    groupname = 'ramses_info'
402
    if(islast) then
403
       groupname = 'ramses_info_last'
404
405
    end if
    call Hdf5_create_group(gr_id,groupname,gr_ramses_id)
406
    aname = 'ncpu'
407
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%ncpu)
408
    aname = 'ndim'
409
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%ndim)
410
    aname = 'levelmin'
411
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%levelmin)
412
    aname = 'levelmax'
413
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%levelmax)
414
    aname = 'ngridmax'
415
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%ngridmax)
416
    aname = 'nstep_coarse'
417
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%nstep_coarse)
418
    aname = 'boxlen'
419
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%boxlen)
420
    aname = 'time'
421
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%time)
422
    aname = 'aexp'
423
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%aexp)
424
    aname = 'h0'
425
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%h0)
426
    aname = 'omega_m'
427
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%omega_m)
428
    aname = 'omega_l'
429
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%omega_l)
430
    aname = 'omega_k'
431
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%omega_k)
432
    aname = 'omega_b'
433
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%omega_b)
434
    aname = 'unit_l'
435
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%unit_l)
436
    aname = 'unit_d'
437
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%unit_d)
438
    aname = 'unit_t'
439
440
441
442
443
444
445
    call Hdf5_write_attr(gr_ramses_id,aname,inforamses%unit_t)
    call Hdf5_close_group(gr_ramses_id)
    call Hdf5_close_group(gr_id)

#ifdef DEBUG
    write(OUTPUT_UNIT,*) 'Exit Write_meta_info_ramses on process ', procid
#endif
446

447
  end subroutine Write_meta_info_ramses
448
449
450
451


  !=======================================================================
  !> Write light cone info as metadata
452
  subroutine Write_meta_info_cone(file_id, infocone, islast)
453

454
455
456
457
458
459
    integer(kind=hid_t), intent(in) :: file_id
    class(type_info_cone), intent(in) :: infocone
    logical(kind=4), intent(in) :: islast

    character(len=H5_STR_LEN) :: groupname
    character(len=H5_STR_LEN) :: aname
460

461
462
    integer(kind=hid_t) :: gr_id
    integer(kind=hid_t) :: gr_cone_id
463

464
465
466
#ifdef DEBUG
    write(OUTPUT_UNIT,*) 'Enter Write_meta_info_cone on process ', procid
#endif
467
468

    groupname = 'metadata'
469
    call Hdf5_open_group(file_id, groupname, gr_id)
470
471

    groupname = 'cone_info'
472
    if(islast) then
473
       groupname = 'cone_info_last'
474
475
    end if
    call Hdf5_create_group(gr_id, groupname, gr_cone_id)
476
    aname = 'ncpu'
477
    call Hdf5_write_attr(gr_cone_id, aname, infocone%ncpu)
478
    aname = 'nstride'
479
    call Hdf5_write_attr(gr_cone_id, aname, infocone%nstride)
480
    aname = 'nstep_coarse'
481
    call Hdf5_write_attr(gr_cone_id, aname, infocone%nstep_coarse)
482
    aname = 'cone_id'
483
    call Hdf5_write_attr(gr_cone_id, aname, infocone%cone_id)
484
    aname = 'nglobalfile'
485
    call Hdf5_write_attr(gr_cone_id, aname, infocone%nglobalfile)
486
    aname = 'isfullsky'
487
    call Hdf5_write_attr(gr_cone_id, aname, infocone%isfullsky)
488
    aname = 'aexp'
489
    call Hdf5_write_attr(gr_cone_id, aname, infocone%aexp)
490
    aname = 'observer_x'
491
    call Hdf5_write_attr(gr_cone_id, aname, infocone%observer_x)
492
    aname = 'observer_y'
493
    call Hdf5_write_attr(gr_cone_id, aname, infocone%observer_y)
494
    aname = 'observer_z'
495
    call Hdf5_write_attr(gr_cone_id, aname, infocone%observer_z)
496
    aname = 'observer_rds'
497
    call Hdf5_write_attr(gr_cone_id, aname, infocone%observer_rds)
498
    aname = 'cone_zlim'
499
    call Hdf5_write_attr(gr_cone_id, aname, infocone%cone_zlim)
500
    aname = 'amax'
501
    call Hdf5_write_attr(gr_cone_id, aname, infocone%amax)
502
    aname = 'amin'
503
    call Hdf5_write_attr(gr_cone_id, aname, infocone%amin)
504
    aname = 'zmax'
505
    call Hdf5_write_attr(gr_cone_id, aname, infocone%zmax)
506
    aname = 'zmin'
507
    call Hdf5_write_attr(gr_cone_id, aname, infocone%zmin)
508
    aname = 'dmax'
509
    call Hdf5_write_attr(gr_cone_id, aname, infocone%dmax)
510
    aname = 'dmin'
511
    call Hdf5_write_attr(gr_cone_id, aname, infocone%dmin)
512
    aname = 'dtol'
513
    call Hdf5_write_attr(gr_cone_id, aname, infocone%dtol)
514
    aname = 'thetay'
515
    call Hdf5_write_attr(gr_cone_id, aname, infocone%thetay)
516
    aname = 'thetaz'
517
    call Hdf5_write_attr(gr_cone_id, aname, infocone%thetaz)
518
    aname = 'theta'
519
    call Hdf5_write_attr(gr_cone_id, aname, infocone%theta)
520
    aname = 'phi'
521
    call Hdf5_write_attr(gr_cone_id, aname, infocone%phi)
522
    aname = 'aendconem2'
523
    call Hdf5_write_attr(gr_cone_id, aname, infocone%aendconem2)
524
    aname = 'aendconem1'
525
    call Hdf5_write_attr(gr_cone_id, aname, infocone%aendconem1)
526
    aname = 'aendcone'
527
    call Hdf5_write_attr(gr_cone_id, aname, infocone%aendcone)
528
    aname = 'zendconem2'
529
    call Hdf5_write_attr(gr_cone_id, aname, infocone%zendconem2)
530
    aname = 'zendconem1'
531
    call Hdf5_write_attr(gr_cone_id, aname, infocone%zendconem1)
532
    aname = 'zendcone'
533
    call Hdf5_write_attr(gr_cone_id, aname, infocone%zendcone)
534
    aname = 'dendconem2'
535
    call Hdf5_write_attr(gr_cone_id, aname, infocone%dendconem2)
536
    aname = 'dendconem1'
537
    call Hdf5_write_attr(gr_cone_id, aname, infocone%dendconem1)
538
    aname = 'dendcone'
539
    call Hdf5_write_attr(gr_cone_id, aname, infocone%dendcone)
540
    aname = 'future'
541
    call Hdf5_write_attr(gr_cone_id, aname, infocone%future)
542
543
544
    select type(infocone)
    type is(type_info_cone_part)
       aname = 'npart'
545
       call Hdf5_write_data(gr_cone_id, aname, infocone%npart)
546
       aname = 'aexpold'
547
       call Hdf5_write_attr(gr_cone_id, aname, infocone%aexpold)
548
       aname = 'zexpold'
549
       call Hdf5_write_attr(gr_cone_id, aname, infocone%zexpold)
550
       aname = 'zexp'
551
       call Hdf5_write_attr(gr_cone_id, aname, infocone%zexp)
552
       aname = 'dexpold'
553
       call Hdf5_write_attr(gr_cone_id, aname, infocone%dexpold)
554
       aname = 'dexp'
555
       call Hdf5_write_attr(gr_cone_id, aname, infocone%dexp)
556
557
    type is(type_info_cone_grav)
       aname = 'nglobalcell'
558
       call Hdf5_write_data(gr_cone_id, aname, infocone%nglobalcell)
559
       aname = 'nlevel'
560
       call Hdf5_write_attr(gr_cone_id, aname, infocone%nlevel)
561
       aname = 'levelmin'
562
       call Hdf5_write_attr(gr_cone_id, aname, infocone%levelmin)
563
       aname = 'levelmax'
564
       call Hdf5_write_attr(gr_cone_id, aname, infocone%levelmax)
565
    end select
566
567
568
569
570
571
    call Hdf5_close_group(gr_cone_id)
    call Hdf5_close_group(gr_id)

#ifdef DEBUG
    write(OUTPUT_UNIT,*) 'Exit Write_meta_info_cone on process ', procid
#endif
572

573
  end subroutine Write_meta_info_cone
574

575
end module modwritemeta