Как получить вывод root .iterate (), идентичный выводу root_ numpy root2array () с неровными массивами - PullRequest
Новая
       8

Как получить вывод root .iterate (), идентичный выводу root_ numpy root2array () с неровными массивами

1 голос
/ 14 марта 2020

Я знаю, что есть решение для аналогичного вопроса: Как получить root .iterate () вывод как root_ numpy root2array () вывод быстро Но, как я понимаю, это подходит только для квартиры ROOT TT. Я хочу иметь обобщенное решение для:

  1. фиксированного размера , но вложенных данных, таких как импульс частицы (px, py, pz), которые представлены в дереве ROOT как vector<double>
  2. данные измерений произвольного размера

Все мои попытки применить asjagged не увенчались успехом. Можно ли избежать jaggedarray для случая (1)?

1 Ответ

1 голос
/ 14 марта 2020

Если данные имеют фиксированный размер, но хранятся как vector<double>, то они обрабатываются так, как если бы они не были фиксированного размера. Up root всегда будет читать их как неровные массивы, и поэтому метод asarray, описанный в другом вопросе, недоступен.

Тем не менее, если у вас больше знаний, чем у метаданных вашего файла, и вы хотите попробовать неподдерживаемый хак, вы можете заставить интерпретацию вашего vector<double> всегда иметь три элемента. Вот пример - во-первых, нам нужен подходящий ROOT файл: 

TFile *f = new TFile("tmp.root", "recreate");
TTree *t = new TTree("t", "t");
std::vector<double> x;
t->Branch("x", &x);
x = {101, 102, 103};
t->Fill();
x = {201, 202, 203};
t->Fill();
x = {301, 302, 303};
t->Fill();
x = {401, 402, 403};
t->Fill();
x = {501, 502, 503};
t->Fill();
x = {601, 602, 603};
t->Fill();
t->Write();
f->Close();

Теперь, обычный способ прочитать это в Up root с помощью JaggedArray: 

>>> import uproot
>>> branch = uproot.open("tmp.root")["t"]["x"]
>>> branch.array()
<JaggedArray [[101.0 102.0 103.0] [201.0 202.0 203.0] [301.0 302.0 303.0] [401.0 402.0 403.0] [501.0 502.0 503.0] [601.0 602.0 603.0]] at 0x7f325ab4eb90>
>>> branch.interpretation
asjagged(asdtype('>f8'), 10)

Но это выделяет новые массивы каждый раз, когда вы читаете его, и выполняет некоторые манипуляции, чтобы найти границы, которые, как вы знаете, являются регулярными. (И до root этого не знает.)

Длина заголовка для вектора STL составляет 10 байтов. После этого его содержимое сериализуется по порядку, от первого до последнего, прежде чем перейти к следующему вектору STL. Для трех 8-байтовых чисел с плавающей запятой это 10 + 8 + 8 + 8 = 34 байта, всегда одинакового размера. Тот факт, что он всегда имеет одинаковый размер, имеет решающее значение для следующего.

A NumPy структурированный массив может представлять массивы структур фиксированного размера в виде dtype: 

>>> array = branch.array(uproot.asdtype([("header", "S10"), ("x", ">f8"), ("y", ">f8"), ("z", ">f8")]))
>>> array
array([(b'@\x00\x00\x1e\x00\t\x00\x00\x00\x03', 101., 102., 103.),
       (b'@\x00\x00\x1e\x00\t\x00\x00\x00\x03', 201., 202., 203.),
       (b'@\x00\x00\x1e\x00\t\x00\x00\x00\x03', 301., 302., 303.),
       (b'@\x00\x00\x1e\x00\t\x00\x00\x00\x03', 401., 402., 403.),
       (b'@\x00\x00\x1e\x00\t\x00\x00\x00\x03', 501., 502., 503.),
       (b'@\x00\x00\x1e\x00\t\x00\x00\x00\x03', 601., 602., 603.)],
      dtype=[('header', 'S10'), ('x', '<f8'), ('y', '<f8'), ('z', '<f8')])

Нам не нравится смотреть на заголовок, поэтому мы можем вырезать его, используя обычную NumPy вырезку: 

>>> array[["x", "y", "z"]]
array([(101., 102., 103.), (201., 202., 203.), (301., 302., 303.),
       (401., 402., 403.), (501., 502., 503.), (601., 602., 603.)],
      dtype={'names':['x','y','z'], 'formats':['<f8','<f8','<f8'], 'offsets':[10,18,26], 'itemsize':34})

(или просто попросить array["x"], чтобы получить не структурированный массив). Поскольку мы можем сделать это с простым uproot.asdtype, мы можем сделать это с предварительно выделенным uproot.asarray. 

>>> import numpy as np
>>> big = np.empty(1000, dtype=[("header", "S10"), ("x", ">f8"), ("y", ">f8"), ("z", ">f8")])
>>> branch.array(uproot.asarray(big.dtype, big))
array([(b'@\x00\x00\x1e\x00\t\x00\x00\x00\x03', 101., 102., 103.),
       (b'@\x00\x00\x1e\x00\t\x00\x00\x00\x03', 201., 202., 203.),
       (b'@\x00\x00\x1e\x00\t\x00\x00\x00\x03', 301., 302., 303.),
       (b'@\x00\x00\x1e\x00\t\x00\x00\x00\x03', 401., 402., 403.),
       (b'@\x00\x00\x1e\x00\t\x00\x00\x00\x03', 501., 502., 503.),
       (b'@\x00\x00\x1e\x00\t\x00\x00\x00\x03', 601., 602., 603.)],
      dtype=[('header', 'S10'), ('x', '>f8'), ('y', '>f8'), ('z', '>f8')])

. Выше, big должен быть по крайней мере настолько большим, насколько набор данных вы хотите чтение и чтение возвращают усеченное представление этого массива. Он не выделяет новый массив - данные хранятся в big: 

>>> big[["x", "y", "z"]][:10]
array([( 1.01000000e+002,  1.02000000e+002,  1.03000000e+002),
       ( 2.01000000e+002,  2.02000000e+002,  2.03000000e+002),
       ( 3.01000000e+002,  3.02000000e+002,  3.03000000e+002),
       ( 4.01000000e+002,  4.02000000e+002,  4.03000000e+002),
       ( 5.01000000e+002,  5.02000000e+002,  5.03000000e+002),
       ( 6.01000000e+002,  6.02000000e+002,  6.03000000e+002),
       ( 1.22164945e-309,  5.26335088e-310,  1.22167067e-309),
       (-5.70498984e+158,  5.97958175e+158, -5.97958175e+158),
       (-4.92033505e+032, -4.92033505e+032, -4.92033505e+032),
       ( 3.77957352e+011,  3.77957221e+011,  3.77957320e+011)],
      dtype={'names':['x','y','z'], 'formats':['>f8','>f8','>f8'], 'offsets':[10,18,26], 'itemsize':34})

Все, кроме 6 прочитанных событий, является неинициализированным мусором, поэтому я рекомендую использовать усеченный вывод из * Функция 1037 *; это просто для того, чтобы показать, что big заполняется - вы не получаете новый массив.

Согласно вашему вопросу (2), если данные не регулярны (в количестве байтов на запись) тогда ты не сможешь это сделать. Кроме того, помните, что вышеупомянутая методика официально не поддерживается: вы должны знать, что ваши данные являются регулярными, и вы должны знать, что векторы STL имеют 10-байтовый заголовок, чего мы не ожидаем от большинства пользователей.

...