openPMD / openPMD-api

:floppy_disk: C++ & Python API for Scientific I/O
https://openpmd-api.readthedocs.io
GNU Lesser General Public License v3.0
134 stars 51 forks source link

Read particle data from openpmd output #1448

Closed StevE-Ong closed 1 year ago

StevE-Ong commented 1 year ago

Hi, I am trying to access the particle data created from picongpu using the following script. But I am not sure it read correctly. The result of the particle position is weird with one very large value and the other very small. The weighting also cannot be read. I appreciate any one could point out the correct reading of particle data. Thanks.

import openpmd_api as io

series = io.Series("./gas_126175.h5",
                io.Access.read_only)

print("Read a Series with openPMD standard version %s" %
        series.openPMD)

i = series.iterations[126175]

for m in i.meshes:
        print("\t {0}".format(m))
        print("")
        print("Iteration X contains {0} particle species:".format(
                len(i.particles)))
for ps in i.particles:
        print("\t {0}".format(ps))
        print("With records:")
        for r in i.particles[ps]:
                print("\t {0}".format(r))

# printing a scalar value
electrons = i.particles["e_all"]
charge = electrons["charge"][io.Mesh_Record_Component.SCALAR]

x = electrons["position"]["x"][()]

print(x)

y = electrons["position"]["y"][()]
z = electrons["position"]["z"][()]
ux = electrons["momentum"]["x"][()] 
uy = electrons["momentum"]["y"][()]
uz = electrons["momentum"]["z"][()]
w = electrons["weighting"]

series.flush()

OUTPUT

Read a Series with openPMD standard version 1.1.0
     picongpu_idProvider

Iteration X contains 1 particle species:
     e_all
With records:
     charge
     mass
     momentum
     momentumPrev1
     position
     positionOffset
     radiationMask
     weighting
[-2.0551398e+08  3.0786527e-41  4.3858892e-15 ... -4.6660318e+00
 -4.0074358e+00 -4.1378455e+00]
Series constructor called with explicit iteration suggests loading a single file with groupBased iteration encoding. Loaded file is fileBased.
---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
[/home/ong/Jupyter/openPMD2dat/openPMD2fortran_2.ipynb](https://file+.vscode-resource.vscode-cdn.net/home/ong/Jupyter/openPMD2dat/openPMD2fortran_2.ipynb) Cell 1 in 
     [33](vscode-notebook-cell:/home/ong/Jupyter/openPMD2dat/openPMD2fortran_2.ipynb#W0sZmlsZQ%3D%3D?line=32) uy = electrons["momentum"]["y"][()]
     [34](vscode-notebook-cell:/home/ong/Jupyter/openPMD2dat/openPMD2fortran_2.ipynb#W0sZmlsZQ%3D%3D?line=33) uz = electrons["momentum"]["z"][()]
---> [35](vscode-notebook-cell:/home/ong/Jupyter/openPMD2dat/openPMD2fortran_2.ipynb#W0sZmlsZQ%3D%3D?line=34) w = electrons["weighting"][()]
     [37](vscode-notebook-cell:/home/ong/Jupyter/openPMD2dat/openPMD2fortran_2.ipynb#W0sZmlsZQ%3D%3D?line=36) series.flush()

TypeError: __getitem__(): incompatible function arguments. The following argument types are supported:
    1. (self: openpmd_api.openpmd_api_cxx.Record_Component_Container, arg0: str) -> openPMD::RecordComponent

Invoked with: , ()
franzpoeschel commented 1 year ago

Hello @StevE-Ong,

I see two small issues:

  1. Data is only read when explicitly flushing. So before printing the data, you would need to insert a series.flush() call, else you get uninitialized memory.

    x = electrons["position"]["x"][()]
    # no data yet in x
    series.flush()
    # now, x is filled
    print(x)
  2. Much like the charge, also the scaling is a scalar record, i.e. you need to write something like w = electrons["weighting"][io.Record_Component.SCALAR][()]. (This is a confusing API design on our side and will hopefully no longer be necessary once #1154 is finished)

Otherwise this looks good at first glance. Do tell if you run into any further problems.

StevE-Ong commented 1 year ago

Hi @franzpoeschel,

Thanks. It works fine now.

Just another small question. Is it possible to read openpmd input using fortran?

franzpoeschel commented 1 year ago

Is it possible to read openpmd input using fortran?

Not with the openPMD-api, no. The best alternatives are:

  1. Use the Fortran bindings of ADIOS2 or HDF5, they have them
  2. Or write a small C++-Fortran interface that provides the functionality which you need