Uploaded image for project: 'ROOT'
  1. ROOT
  2. ROOT-10148

[TTreReader] Cannot read columns holding arrays of opaque types (Double32_t and Float16_t)

    XMLWordPrintable

    Details

    • Type: Bug
    • Status: Closed (View Workflow)
    • Priority: High
    • Resolution: Completed
    • Affects Version/s: 6.16/00
    • Fix Version/s: 6.18/00
    • Component/s: I/O
    • Labels:
      None
    • Environment:

      all

      Description

      This program* illustrates the problem. It can be compiled and run following the instructions at its beginning.
      One can read the columns using the bare TTree interface but TTreeReader is not able to do so, causing RDF not to be able too either.

      *

      // Compilation and execution:
      // To test float16:
      // g++ -D FLOAT16T -o testfloat16 testOpaqueTypesNanoAOD.cpp `root-config
      // --cflags --libs`
      // ./testfloat16 or ./testfloat16 "arr[n]/[0,2pi,8]" or ./testfloat16 "arr[n]/f[0,twopi,10]"
      // To test with float:
      // g++ -o testfloat testOpaqueTypesNanoAOD.cpp `root-config --cflags --libs`
      // ./testfloat
       
      #include <ROOT/TSeq.hxx>
      #include <TFile.h>
      #include <TMath.h>
      #include <TRandom.h>
      #include <TTree.h>
       
      #include <iostream>
       
      const auto nevts = 100000;
      const auto treename = "tree";
       
      #ifdef FLOAT16T
      using flType = Float16_t;
      const auto leaflist = "arr[n]/f";
      const auto filename = "out_float16.root";
      #else
      using flType = float;
      const auto leaflist = "arr[n]/F";
      const auto filename = "out_float.root";
      #endif
       
      void write(const char *theLeafList) {
        std::cout << "Using leaflist \"" << theLeafList << "\"" << std::endl;
        TFile f(filename, "RECREATE");
        TTree t("t", "t");
        int n;
        flType arr[64];
        t.Branch("n", &n);
        t.Branch("arr", arr, theLeafList);
        for (auto ievt : ROOT::TSeqU(nevts)) {
          n = gRandom->Uniform(0, 32);
          for (auto inumb : ROOT::TSeqI(n)) {
            auto rndm = gRandom->Uniform(0, 2 * TMath::Pi());
            arr[inumb] = rndm;
          }
          t.Fill();
        }
        t.Write();
      }
       
      void read() {
        // old school reading
        flType arr[64];
        int n;
        TFile f(filename);
        auto t = f.Get<TTree>("t");
        t->SetBranchAddress("arr", arr);
        t->SetBranchAddress("n", &n);
        for (auto ievt : ROOT::TSeqU(100000)) {
          t->GetEntry();
          for (auto i : ROOT::TSeqU(n)) {
            std::cout << "Evt " << ievt << " v = " << arr[i] << std::endl;
          }
        }
      }
       
      int main(int argc, char **argv) {
        gRandom->SetSeed(1);
        const char *theLeafList;
      #ifdef FLOAT16T
        if (argc == 1) {
          theLeafList = leaflist;
        } else {
          theLeafList = argv[1];
        }
      #else
        theLeafList = leaflist;
      #endif
       
        write(theLeafList);
        read();
        return 0;
      }
      
      

        Attachments

          Activity

            People

            • Assignee:
              dpiparo Danilo Piparo
              Reporter:
              dpiparo Danilo Piparo
            • Votes:
              0 Vote for this issue
              Watchers:
              1 Start watching this issue

              Dates

              • Created:
                Updated:
                Resolved:
                Actual End: