go-hep.org/x/hep@v0.38.1/groot/riofs/gendata/gen-ndim-slice.go (about)

     1  // Copyright ©2020 The go-hep Authors. All rights reserved.
     2  // Use of this source code is governed by a BSD-style
     3  // license that can be found in the LICENSE file.
     4  
     5  //go:build ignore
     6  
     7  package main
     8  
     9  import (
    10  	"flag"
    11  	"log"
    12  
    13  	"go-hep.org/x/hep/groot/internal/rtests"
    14  )
    15  
    16  var (
    17  	root  = flag.String("f", "test-ndim-slice.root", "output ROOT file")
    18  	split = flag.Int("split", 0, "default split-level for TTree")
    19  )
    20  
    21  func main() {
    22  	flag.Parse()
    23  
    24  	out, err := rtests.RunCxxROOT("gentree", []byte(script), *root, *split)
    25  	if err != nil {
    26  		log.Fatalf("could not run ROOT macro:\noutput:\n%v\nerror: %+v", string(out), err)
    27  	}
    28  }
    29  
    30  const script = `
    31  #include <string.h>
    32  #include <stdio.h>
    33  
    34  #define MAXSLICE 5
    35  
    36  #define OFFSET 0
    37  #define N1 2
    38  #define N2 3
    39  #define N3 4
    40  
    41  struct Event {
    42  	int32_t  N;
    43  
    44  	bool     SliceBs[MAXSLICE][2][3][4];  //[N]
    45  	int8_t   SliceI8[MAXSLICE][2][3][4];  //[N]
    46  	int16_t  SliceI16[MAXSLICE][2][3][4]; //[N]
    47  	int32_t  SliceI32[MAXSLICE][2][3][4]; //[N]
    48  	int64_t  SliceI64[MAXSLICE][2][3][4]; //[N]
    49  	uint8_t  SliceU8[MAXSLICE][2][3][4];  //[N]
    50  	uint16_t SliceU16[MAXSLICE][2][3][4]; //[N]
    51  	uint32_t SliceU32[MAXSLICE][2][3][4]; //[N]
    52  	uint64_t SliceU64[MAXSLICE][2][3][4]; //[N]
    53  	float    SliceF32[MAXSLICE][2][3][4]; //[N]
    54  	double   SliceF64[MAXSLICE][2][3][4]; //[N]
    55  
    56  	Float16_t    SliceD16[MAXSLICE][2][3][4]; //[N]
    57  	Double32_t   SliceD32[MAXSLICE][2][3][4]; //[N]
    58  
    59  };
    60  
    61  void gentree(const char* fname, int splitlvl = 99) {
    62  	int bufsize = 32000;
    63  	int evtmax = 2;
    64  
    65  	auto f = TFile::Open(fname, "RECREATE");
    66  	auto t = new TTree("tree", "my tree title");
    67  
    68  	Event e;
    69  
    70  	t->Branch("N", &e.N, "N/I");
    71  
    72  	t->Branch("SliBs",  e.SliceBs,  "SliBs[N][2][3][4]/O");
    73  	t->Branch("SliI8",  e.SliceI8,  "SliI8[N][2][3][4]/B");
    74  	t->Branch("SliI16", e.SliceI16, "SliI16[N][2][3][4]/S");
    75  	t->Branch("SliI32", e.SliceI32, "SliI32[N][2][3][4]/I");
    76  	t->Branch("SliI64", e.SliceI64, "SliI64[N][2][3][4]/L");
    77  	t->Branch("SliU8",  e.SliceU8,  "SliU8[N][2][3][4]/b");
    78  	t->Branch("SliU16", e.SliceU16, "SliU16[N][2][3][4]/s");
    79  	t->Branch("SliU32", e.SliceU32, "SliU32[N][2][3][4]/i");
    80  	t->Branch("SliU64", e.SliceU64, "SliU64[N][2][3][4]/l");
    81  	t->Branch("SliF32", e.SliceF32, "SliF32[N][2][3][4]/F");
    82  	t->Branch("SliF64", e.SliceF64, "SliF64[N][2][3][4]/D");
    83  
    84  	t->Branch("SliD16", e.SliceD16, "SliD16[N][2][3][4]/f[0,0,16]");
    85  	t->Branch("SliD32", e.SliceD32, "SliD32[N][2][3][4]/d[0,0,32]");
    86  
    87  	for (int j = 0; j != evtmax; j++) {
    88  		int i = j + OFFSET;
    89  		e.N = i+1;
    90  		for (int i0 = 0; i0 != e.N; i0++) {
    91  			for (int i1 = 0; i1 != N1; i1++) {
    92  				for (int i2 = 0; i2 != N2; i2++) {
    93  					for (int i3 = 0; i3 != N3; i3++) {
    94  						e.SliceBs[i0][i1][i2][i3]  = i3%2 == 0;
    95  						e.SliceI8[i0][i1][i2][i3]  = -i;
    96  						e.SliceI16[i0][i1][i2][i3] = -i;
    97  						e.SliceI32[i0][i1][i2][i3] = -i;
    98  						e.SliceI64[i0][i1][i2][i3] = -i;
    99  						e.SliceU8[i0][i1][i2][i3]  = i;
   100  						e.SliceU16[i0][i1][i2][i3] = i;
   101  						e.SliceU32[i0][i1][i2][i3] = i;
   102  						e.SliceU64[i0][i1][i2][i3] = i;
   103  						e.SliceF32[i0][i1][i2][i3] = float(i);
   104  						e.SliceF64[i0][i1][i2][i3] = double(i);
   105  						e.SliceD16[i0][i1][i2][i3] = float(i);
   106  						e.SliceD32[i0][i1][i2][i3] = double(i);
   107  						i++;
   108  					}
   109  				}
   110  			}
   111  		}
   112  		t->Fill();
   113  	}
   114  
   115  	f->Write();
   116  	f->Close();
   117  
   118  	exit(0);
   119  }
   120  `