Skip to content

Quick Start

Reading a HIPO File

The most common workflow: open a file, load its schema dictionary, create banks, and iterate through events.

#include "reader.h"

int main() {
    // 1. Open the file and load schemas
    hipo::reader reader;
    reader.open("data.hipo");

    hipo::dictionary dict;
    reader.readDictionary(dict);

    // 2. Create bank objects from schemas
    hipo::bank particles(dict.getSchema("REC::Particle"));
    hipo::event event;

    // 3. Iterate through events
    while (reader.next()) {
        reader.read(event);
        event.read(particles);

        for (int row = 0; row < particles.getRows(); row++) {
            int   pid = particles.getInt("pid", row);
            float px  = particles.getFloat("px", row);
            float py  = particles.getFloat("py", row);
            float pz  = particles.getFloat("pz", row);
            printf("%d: pid=%d  p=(%.3f, %.3f, %.3f)\n",
                   row, pid, px, py, pz);
        }
    }
}

Banklist Shortcut

For multiple banks, use the banklist convenience API:

hipo::banklist list = reader.getBanks({"REC::Particle", "REC::Event"});

while (reader.next(list)) {
    // list[0] is REC::Particle, list[1] is REC::Event
    int nrows = list[0].getRows();
    // ...
}

Writing a HIPO File

#include "writer.h"

int main() {
    // 1. Define schemas
    hipo::schema schemaPart("event::particle", 100, 1);
    schemaPart.parse("pid/S,px/F,py/F,pz/F");

    // 2. Create writer and register schemas BEFORE opening
    hipo::writer writer;
    writer.getDictionary().addSchema(schemaPart);
    writer.open("output.hipo");

    // 3. Create events with banks
    hipo::event event;
    for (int i = 0; i < 1000; i++) {
        int nrows = 2 + rand() % 10;
        hipo::bank bank(schemaPart, nrows);

        for (int row = 0; row < nrows; row++) {
            bank.putShort("pid", row, 211);
            bank.putFloat("px", row, 0.5f * row);
            bank.putFloat("py", row, 0.3f * row);
            bank.putFloat("pz", row, 1.0f + 0.1f * row);
        }

        event.reset();
        event.addStructure(bank);
        writer.addEvent(event);
    }

    writer.close();
}

Important

Schemas must be added to the writer's dictionary before calling open(). They are written into the file header.

Type System

HIPO supports 6 data types, specified by single characters in schema definitions:

Character Type C++ Type Size
B Byte int8_t 1 byte
S Short int16_t 2 bytes
I Int int32_t 4 bytes
F Float float 4 bytes
D Double double 8 bytes
L Long int64_t 8 bytes

Schema format: "column1/T,column2/T,..." where T is the type character.

hipo::schema s("REC::Particle", 300, 1);
s.parse("pid/I,px/F,py/F,pz/F,vx/D,vy/D,vz/D,charge/B,status/L");

Python Quick Start

The hipopy package mirrors the C++ API with Pythonic conventions and NumPy integration.

Reading a HIPO File

import hipopy
import numpy as np

# High-level API with context manager
with hipopy.open("data.hipo", banks=["REC::Particle"]) as f:
    for banks in f:
        bank = banks["REC::Particle"]
        pid = bank["pid"]       # NumPy int32 array
        px  = bank["px"]        # NumPy float32 array
        py  = bank["py"]
        pt  = np.sqrt(px**2 + py**2)
        print(f"{len(pid)} particles, max pT = {pt.max():.3f}")

The low-level API matches the C++ pattern directly:

reader = hipopy.Reader("data.hipo")
dict_  = hipopy.Dictionary()
reader.read_dictionary(dict_)

bank  = hipopy.Bank(dict_["REC::Particle"])
event = hipopy.Event()

while reader.next():
    reader.read(event)
    event.read(bank)
    for row in range(bank.get_rows()):
        pid = bank.get_int("pid", row)
        px  = bank.get_float("px", row)

Writing a HIPO File

schema = hipopy.Schema("event::particle", 100, 1)
schema.parse("pid/S,px/F,py/F,pz/F")

d = hipopy.Dictionary()
d.add_schema(schema)

writer = hipopy.Writer()
writer.add_dictionary(d)
writer.open("output.hipo")

for i in range(1000):
    event = hipopy.Event()
    bank = hipopy.Bank(schema, 5)
    for row in range(5):
        bank.put_short("pid", row, 211)
        bank.put_float("px", row, 0.5 * row)
        bank.put_float("py", row, 0.3 * row)
        bank.put_float("pz", row, 1.0 + 0.1 * row)
    event.add_structure(bank)
    writer.add_event(event)

writer.close()

Batch Reading

# Flat NumPy arrays (all rows concatenated across events)
data = hipopy.read_columns("data.hipo", "REC::Particle", ["pid", "px", "py"])
pt = np.sqrt(data["px"]**2 + data["py"]**2)

# Ragged Awkward Arrays (per-event structure preserved)
particles = hipopy.read_bank("data.hipo", "REC::Particle", ["pid", "px", "py"])
# particles.px[0] = first event's px values, particles.px[1] = second event's, etc.

See Python API Reference for the complete API.

Next Steps