-
Notifications
You must be signed in to change notification settings - Fork 0
/
ReadEagleSubfind_halo.py
177 lines (159 loc) · 11.6 KB
/
ReadEagleSubfind_halo.py
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
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
import h5py as h5
import os
from numpy import *
def ReadEagleSubfind_halo(Base,DirList,fend,exts):
fn = Base+DirList+'/groups_'+exts+fend+'/'+'eagle_subfind_tab_'+exts+fend+'.0.hdf5'
print ' __' ; print ' Directory:',Base+DirList+'/groups_'+exts+fend+'/'+'eagle_subfind_tab_'+exts+fend+' ...'
fs = h5.File(fn,"r")
Header = fs['Header'].attrs
Ntask = Header['NTask']
TotNgroups = Header['TotNgroups']
TotNsubgroups = Header['TotNsubgroups']
fs.close()
# --- Define Group Arrays
if TotNgroups > 0:
Group_M_Crit200 = empty(TotNgroups, dtype=float)
Group_R_Crit200 = empty(TotNgroups, dtype=float)
#Group_R_Mean200 = empty(TotNgroups, dtype=float)
#GroupMass = empty(TotNgroups, dtype=float)
GroupPos = empty((TotNgroups,3), dtype=float)
#GroupLen = empty(TotNgroups, dtype=int)
FirstSub = empty(TotNgroups, dtype=int)
NumOfSubhalos = empty(TotNgroups, dtype=int)
# --- Define Subhalo Arrays
#HalfMassRad = empty((TotNsubgroups,6), dtype=float)
#HalfMassProjRad = empty((TotNsubgroups,6), dtype=float)
MassType = empty((TotNsubgroups,6), dtype=float)
Pcom = empty((TotNsubgroups,3), dtype=float)
SubPos = empty((TotNsubgroups,3), dtype=float)
#SubLenType = empty((TotNsubgroups,6), dtype=int)
Vbulk = empty((TotNsubgroups,3), dtype=float)
Vmax = empty(TotNsubgroups, dtype=float)
VmaxRadius = empty(TotNsubgroups, dtype=float)
SFRate = empty(TotNsubgroups, dtype=float)
#BHMass = empty(TotNsubgroups, dtype=float)
#StellarMass = empty(TotNsubgroups, dtype=float)
SubMass = empty(TotNsubgroups, dtype=float)
SubLen = empty(TotNsubgroups, dtype=int)
#SubOffset = empty(TotNsubgroups, dtype=int)
SubGroupNumber = empty(TotNsubgroups, dtype=int)
GroupNumber = empty(TotNsubgroups, dtype=int)
#MostBoundID = empty(TotNsubgroups, dtype=int)
Etot = empty(TotNsubgroups, dtype=int)
Ekin = empty(TotNsubgroups, dtype=int)
# Mass_001kpc = empty((TotNsubgroups,6), dtype=float)
# Mass_003kpc = empty((TotNsubgroups,6), dtype=float)
# Mass_005kpc = empty((TotNsubgroups,6), dtype=float)
# Mass_010kpc = empty((TotNsubgroups,6), dtype=float)
# Mass_020kpc = empty((TotNsubgroups,6), dtype=float)
Mass_030kpc = empty((TotNsubgroups,6), dtype=float)
# Mass_040kpc = empty((TotNsubgroups,6), dtype=float)
#Mass_050kpc = empty((TotNsubgroups,6), dtype=float)
# Mass_070kpc = empty((TotNsubgroups,6), dtype=float)
# Mass_100kpc = empty((TotNsubgroups,6), dtype=float)
# SFR_001kpc = empty(TotNsubgroups, dtype=float)
# SFR_003kpc = empty(TotNsubgroups, dtype=float)
# SFR_005kpc = empty(TotNsubgroups, dtype=float)
# SFR_010kpc = empty(TotNsubgroups, dtype=float)
# SFR_020kpc = empty(TotNsubgroups, dtype=float)
SFR_030kpc = empty(TotNsubgroups, dtype=float)
# SFR_040kpc = empty(TotNsubgroups, dtype=float)
# SFR_050kpc = empty(TotNsubgroups, dtype=float)
SFR_070kpc = empty(TotNsubgroups, dtype=float)
# SFR_100kpc = empty(TotNsubgroups, dtype=float)
#subtab_index = empty(TotNsubgroups, dtype=int)
#sub_location = empty(TotNsubgroups, dtype=int)
NGrp_c = 0
NSub_c = 0
for ifile in range(0,Ntask):
if os.path.isdir(Base+DirList+'data/'):
fn = Base+DirList+'data/groups_'+exts+fend+'/'+'eagle_subfind_tab_'+exts+fend+'.'+str(ifile)+'.hdf5'
else:
fn = Base+DirList+'/groups_'+exts+fend+'/'+'eagle_subfind_tab_'+exts+fend+'.'+str(ifile)+'.hdf5'
fs = h5.File(fn,"r")
Header = fs['Header'].attrs
HubbleParam = Header['HubbleParam']
#PartMassDM = Header['PartMassDM']
Redshift = Header['Redshift']
Om = [Header['Omega0'],Header['OmegaLambda'],Header['OmegaBaryon']]
Ngroups = Header['Ngroups']
BoxSize = Header['BoxSize']
#TotNgroups = Header['TotNgroups']
Nsubgroups = Header['Nsubgroups']
#TotNsubgroups = Header['TotNsubgroups']
if Ngroups > 0:
Group_M_Crit200[NGrp_c:NGrp_c+Ngroups] = fs["FOF/Group_M_Crit200"].value
Group_R_Crit200[NGrp_c:NGrp_c+Ngroups] = fs["FOF/Group_R_Crit200"].value
#Group_R_Mean200[NGrp_c:NGrp_c+Ngroups] = fs["FOF/Group_R_Mean200"].value
#GroupLen[NGrp_c:NGrp_c+Ngroups] = fs["FOF/GroupLength"].value
#GroupMass[NGrp_c:NGrp_c+Ngroups] = fs["FOF/GroupMass"].value
GroupPos[NGrp_c:NGrp_c+Ngroups] = fs["FOF/GroupCentreOfPotential"].value
FirstSub[NGrp_c:NGrp_c+Ngroups] = fs["FOF/FirstSubhaloID"].value
NumOfSubhalos[NGrp_c:NGrp_c+Ngroups] = fs["FOF/NumOfSubhalos"].value
NGrp_c += Ngroups
if Nsubgroups > 0:
#HalfMassRad[NSub_c:NSub_c+Nsubgroups,:] = fs["Subhalo/HalfMassRad"].value
#HalfMassProjRad[NSub_c:NSub_c+Nsubgroups,:]= fs["Subhalo/HalfMassProjRad"].value
#SubLenType[NSub_c:NSub_c+Nsubgroups,:] = fs["Subhalo/SubLengthType"].value
MassType[NSub_c:NSub_c+Nsubgroups,:] = fs["Subhalo/MassType"].value
SubPos[NSub_c:NSub_c+Nsubgroups,:] = fs["Subhalo/CentreOfPotential"].value
Pcom[NSub_c:NSub_c+Nsubgroups,:] = fs["Subhalo/CentreOfMass"].value
Vbulk[NSub_c:NSub_c+Nsubgroups] = fs["Subhalo/Velocity"].value
Vmax[NSub_c:NSub_c+Nsubgroups] = fs["Subhalo/Vmax"].value
VmaxRadius[NSub_c:NSub_c+Nsubgroups] = fs["Subhalo/VmaxRadius"].value
SubLen[NSub_c:NSub_c+Nsubgroups] = fs["Subhalo/SubLength"].value
SFRate[NSub_c:NSub_c+Nsubgroups] = fs["Subhalo/StarFormationRate"].value
#BHMass[NSub_c:NSub_c+Nsubgroups] = fs["Subhalo/BlackHoleMass"].value
#StellarMass[NSub_c:NSub_c+Nsubgroups] = fs["Subhalo/Stars/Mass"].value
SubMass[NSub_c:NSub_c+Nsubgroups] = fs["Subhalo/Mass"].value
#SubOffset[NSub_c:NSub_c+Nsubgroups] = fs["Subhalo/SubOffset"].value
SubGroupNumber[NSub_c:NSub_c+Nsubgroups] = fs["Subhalo/SubGroupNumber"].value
GroupNumber[NSub_c:NSub_c+Nsubgroups] = fs["Subhalo/GroupNumber"].value
#MostBoundID[NSub_c:NSub_c+Nsubgroups] = fs["Subhalo/IDMostBound"].value
# Mass_001kpc[NSub_c:NSub_c+Nsubgroups,:] = fs["Subhalo/ApertureMeasurements/Mass/001kpc"].value
# Mass_003kpc[NSub_c:NSub_c+Nsubgroups,:] = fs["Subhalo/ApertureMeasurements/Mass/003kpc"].value
# Mass_005kpc[NSub_c:NSub_c+Nsubgroups,:] = fs["Subhalo/ApertureMeasurements/Mass/005kpc"].value
# Mass_010kpc[NSub_c:NSub_c+Nsubgroups,:] = fs["Subhalo/ApertureMeasurements/Mass/010kpc"].value
# Mass_020kpc[NSub_c:NSub_c+Nsubgroups,:] = fs["Subhalo/ApertureMeasurements/Mass/020kpc"].value
Mass_030kpc[NSub_c:NSub_c+Nsubgroups,:] = fs["Subhalo/ApertureMeasurements/Mass/030kpc"].value
# Mass_040kpc[NSub_c:NSub_c+Nsubgroups,:] = fs["Subhalo/ApertureMeasurements/Mass/040kpc"].value
#Mass_050kpc[NSub_c:NSub_c+Nsubgroups,:] = fs["Subhalo/ApertureMeasurements/Mass/050kpc"].value
# Mass_070kpc[NSub_c:NSub_c+Nsubgroups,:] = fs["Subhalo/ApertureMeasurements/Mass/070kpc"].value
# Mass_100kpc[NSub_c:NSub_c+Nsubgroups,:] = fs["Subhalo/ApertureMeasurements/Mass/100kpc"].value
# SFR_001kpc[NSub_c:NSub_c+Nsubgroups] = fs["Subhalo/ApertureMeasurements/SFR/001kpc"].value
# SFR_003kpc[NSub_c:NSub_c+Nsubgroups] = fs["Subhalo/ApertureMeasurements/SFR/003kpc"].value
# SFR_005kpc[NSub_c:NSub_c+Nsubgroups] = fs["Subhalo/ApertureMeasurements/SFR/005kpc"].value
# SFR_010kpc[NSub_c:NSub_c+Nsubgroups] = fs["Subhalo/ApertureMeasurements/SFR/010kpc"].value
# SFR_020kpc[NSub_c:NSub_c+Nsubgroups] = fs["Subhalo/ApertureMeasurements/SFR/020kpc"].value
SFR_030kpc[NSub_c:NSub_c+Nsubgroups] = fs["Subhalo/ApertureMeasurements/SFR/030kpc"].value
# SFR_040kpc[NSub_c:NSub_c+Nsubgroups] = fs["Subhalo/ApertureMeasurements/SFR/040kpc"].value
# SFR_050kpc[NSub_c:NSub_c+Nsubgroups] = fs["Subhalo/ApertureMeasurements/SFR/050kpc"].value
SFR_070kpc[NSub_c:NSub_c+Nsubgroups] = fs["Subhalo/ApertureMeasurements/SFR/070kpc"].value
# SFR_100kpc[NSub_c:NSub_c+Nsubgroups] = fs["Subhalo/ApertureMeasurements/SFR/100kpc"].value
Etot[NSub_c:NSub_c+Nsubgroups] = fs["Subhalo/TotalEnergy"].value
Ekin[NSub_c:NSub_c+Nsubgroups] = fs["Subhalo/KineticEnergy"].value
#subtab_index[NSub_c:NSub_c+Nsubgroups].fill(ifile)
#sub_location[NSub_c:NSub_c+Nsubgroups] = arange(0,Nsubgroups)
NSub_c += Nsubgroups
fs.close()
return {'Group_M_Crit200':Group_M_Crit200, 'Group_R_Crit200':Group_R_Crit200, 'GroupPos':GroupPos,
#'Group_R_Mean200':Group_R_Mean200, 'GroupMass':GroupMass, #'GroupLen':GroupLen,
'FirstSub':FirstSub, #'subtab_index':subtab_index, 'sub_location':sub_location,
'NumOfSubhalos':NumOfSubhalos, 'Vmax':Vmax, 'Vbulk':Vbulk, 'Rmax':VmaxRadius,
'SubPos':SubPos, 'Pcom':Pcom, 'SubLen':SubLen,
'SubGroupNumber':SubGroupNumber, #'SubOffset':SubOffset, 'HalfMassProjRad':HalfMassProjRad,
'GroupNumber':GroupNumber, #'Mass_001kpc':Mass_001kpc, 'Mass_003kpc':Mass_001kpc,
# 'Mass_003kpc':Mass_003kpc, 'Mass_005kpc':Mass_005kpc, 'Mass_010kpc':Mass_010kpc,
#'Mass_020kpc':Mass_020kpc, 'Mass_070kpc':Mass_070kpc, 'Mass_040kpc':Mass_040kpc,
'Mass_030kpc':Mass_030kpc, 'SFR_030kpc':SFR_030kpc, #'Mass_050kpc':Mass_050kpc,
# 'SFR_001kpc':SFR_001kpc, 'SFR_003kpc':SFR_001kpc, 'Mass_100kpc':Mass_100kpc,
# 'SFR_003kpc':SFR_003kpc, 'SFR_005kpc':SFR_005kpc, 'SFR_010kpc':SFR_010kpc,
'SFR_030kpc':SFR_030kpc, 'SFR_070kpc':SFR_070kpc, #'SFR_040kpc':SFR_040kpc,
# 'SFR_050kpc':SFR_050kpc, 'SFR_020kpc':SFR_020kpc, 'SFR_100kpc':SFR_100kpc,
'TotNgroups':TotNgroups, 'TotNsubgroups':TotNsubgroups,
'Redshift':Redshift, 'BoxSize':BoxSize, 'HubbleParam':HubbleParam,#'PartMassDM':PartMassDM, # 'StellarMass':StellarMass,
'Omega':Om, 'SFRate':SFRate, 'MassType':MassType,
'SubMass':SubMass, 'Etot':Etot, 'Ekin':Ekin # 'MostBoundID':MostBoundID
}
else:
return {'Omega':Om, 'BoxSize':BoxSize, 'Redshift':Redshift, 'HubbleParam':HubbleParam}# 'PartMassDM':PartMassDM}