1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
|
/* Reads in a cube, writes out as root */
/* (c) P.M. Jones */
/* (09/06/20) */
/* Ver. 1.0 */
/* (03/08/20) */
#include <stdio.h>
#include <ctype.h>
#include <string.h>
#include "TFile.h"
#include "TRandom.h"
#include "TH2S.h"
#include "TH3D.h"
#define SIZE 4096
//uint16_t cube[SIZE][SIZE][SIZE]={0};
uint16_t matrix[SIZE][SIZE]={0};
uint16_t matrix2[SIZE][SIZE]={0};
uint16_t xproj[SIZE]={0};
uint16_t yproj[SIZE]={0};
uint16_t zproj[SIZE]={0};
int makecube(int argc, char *argv[])
{ FILE *f;
int i, j, k, k2, z;
char in_name[40];
char out_name[40];
char file_in[150];
if (argc < 2 || argc > 2)
{
fprintf(stderr, "Usage: %s cube\n", argv[0]);
exit(1);
}
strcpy(in_name, argv[1]);
strcpy(out_name, in_name);
strcat(out_name, ".root");
TFile *g = new TFile(out_name,"recreate");
TH1I *xp = new TH1I("xp","Projection X",4096,0,4095);
TH1I *yp = new TH1I("yp","Projection Y",4096,0,4095);
TH1I *zp = new TH1I("zp","Projection Z",4096,0,4095);
TH2S *mat = new TH2S("mat","Matrix",4096,0,4095,4096,0,4095);
TH3D *cub = new TH3D("ggdt_cube", "ggdt cube", 300,0,1500,300,0,1500,600,-3000,3000);
// TH3S *cub = new TH3S("cub","Cube",4096,0,4095,4096,0,4095,4096,0,4095);
/* Read in matrices */
for (k=0; k< 4096; k++)
{
z=k+1;
sprintf(file_in, "/space/lumkile/01Aug2020/%s.%d.mat", in_name, z);
f = fopen(file_in,"rb"); /* open file for binary read */
if (!f)
{
fprintf(stderr, "Can't open file '%s'\n",in_name);
exit(1);
}
printf("Reading matrix '%s' in...\n", file_in);
// fread(cube, sizeof(cube[0][0][k]), SIZE*SIZE, f);
fread(matrix, sizeof(matrix[0][0]), SIZE*SIZE, f);
fclose(f);
for (j=0; j< 3000; j++)
{
for (i=0; i< 3000; i++)
{
matrix2[i][j]=0;
}
}
for (j=0; j< 3000; j++)
{
for (i=0; i< 3000; i++)
{
matrix2[i/2][j/2]=matrix2[i/2][j/2]+matrix[i][j];
xproj[i/2]=xproj[i/2]+matrix[i][j];
yproj[i/2]=yproj[i/2]+matrix[j][i];
mat->Fill(i/2,j/2,matrix2[i/2][j/2]);
k2 = (k-2048);
cub->Fill(i/2,j/2,k2,matrix2[i/2][j/2]);
}
}
}
for (i=0; i< SIZE; i++)
{
xp->Fill(i,xproj[i]);
yp->Fill(i,yproj[i]);
zp->Fill(i,zproj[i]);
}
printf("\nWriting cube '%s' out...", out_name);
printf("\n");
xp->Write();
yp->Write();
zp->Write();
mat->Write();
cub->Write();
delete g;
}
|
