forked from Aramaan/plasma-modelling-using-GTC-X
-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathreadEFITgfile.m
110 lines (88 loc) · 3.57 KB
/
readEFITgfile.m
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
function A = readEFITgfile(gfile)
% open file for reading
findex = fopen(gfile,'r');
if findex < 0 %error
error(['Could not open the gfile ', gfile])
end
% read data
% header, pest parameter, # horizontal R grid points, # vertical z grid points
% READ(kueql,1000)(etitl(i),i=1,5),date,ipestg,nx,nz
A.case = fscanf(findex,'%48c',1);
A.idum = fscanf(findex,'%d',1);
A.nw = fscanf(findex,'%d',1);
A.nh = fscanf(findex,'%d',1);
% horz size (m) of computational box, vert size (m) of computational box, R of B location, R of boundary, middle z value
% READ(kueql,1010) xdim,zdim,rcnt,redge,zmid
A.rdim = fscanf(findex, '%e', 1);
A.zdim = fscanf(findex, '%e', 1);
A.rcentr = fscanf(findex, '%e', 1);
A.rleft = fscanf(findex, '%e', 1);
A.zmid = fscanf(findex, '%e', 1);
% major radius (at magnetic axis), z at mag axis, polflux at mag axis, polflux at limiter, B at specified R location (see above)
% READ(kueql,1010) xma,zma,psimax,psilim,btor
A.rmaxis = fscanf(findex, '%e', 1);
A.zmaxis = fscanf(findex, '%e', 1);
A.simag = fscanf(findex, '%e', 1);
A.sibry = fscanf(findex, '%e', 1);
A.bcentr = fscanf(findex, '%e', 1);
% current in Amps,
% READ(kueql,1010) totcur,psimx(1),psimx(2),xax(1),xax(2)
A.current = fscanf(findex, '%e', 1);
tmp = fscanf(findex, '%e', 2);
tmp = fscanf(findex, '%e', 2);
% READ(kueql,1010) zax(1),zax(2),psisep,xsep,zsep
tmp = fscanf(findex, '%e', 2);
tmp = fscanf(findex, '%e', 1);
tmp = fscanf(findex, '%e', 1);
A.xdum = fscanf(findex, '%e', 1);
%poloidal current function... R*B_toroidal (Wb/m)
% READ(kueql,1010) (sf(i), i = 1,npv)
A.fpol = fscanf(findex, '%e', A.nw);
%plasma pressure in Joule/m^3
% READ(kueql,1010) (sp(i), i = 1,npv)
A.pres = fscanf(findex, '%e', A.nw);
%f df/dpsi
% READ(kueql,1010) (sffp(i), i = 1,npv)
A.ffprim = fscanf(findex, '%e', A.nw);
% dp/dpsi
% READ(kueql,1010) (spp(i), i = 1,npv)
A.pprime = fscanf(findex, '%e', A.nw);
% polflux on rectangular grid
% READ(kueql,1010) ((psi(i,j), i = 1,nx), j = 1,nz)
tmppsi = fscanf(findex, '%e', A.nw*A.nh);
A.psizr = reshape(tmppsi,A.nw,A.nh);
% q values on flux surfaces
% READ(kueql,1010) (efitq(i), i = 1,npv)
A.qpsi = fscanf(findex, '%e', A.nw);
% number of boundary points, and number of limiter points
% READ(kueql,1020) neftbd,neftlm
A.nbbbs = fscanf(findex, '%d', 1);
A.limitr = fscanf(findex, '%d', 1);
% plot the boundary
% READ(kueql,1010) (efitbdy(i),i = 1,nefbdy)
tmpbdy = fscanf(findex, '%e', 2*A.nbbbs);
tmpbdy = reshape(tmpbdy,2,A.nbbbs);
A.rbbbs = tmpbdy(1,:);
A.zbbbs = tmpbdy(2,:);
% plot the limiter
% READ(kueql,1010) (efitlim(i),i = 1,neflim)
tmplim = fscanf(findex, '%e', 2*A.limitr);
tmplim = reshape(tmplim,2,A.limitr);
A.rlim = tmplim(1,:);
A.zlim = tmplim(2,:);
% check some switches, and characteristic radius of rotation
A.kvtor = fscanf(findex, '%d', 1);
A.rvtor = fscanf(findex, '%e', 1);
A.nmass = fscanf(findex, '%d', 1);
% rotation related quantities; optional
if A.kvtor > 0
A.pressw = fscanf(findex, '%e', A.nw);
A.pwprim = fscanf(findex, '%e', A.nw);
end
% mass related properties; optional
if A.nmass > 0
A.dmion = fscanf(findex, '%e', A.nw);
end
% toroidal flux
A.rhovn = fscanf(findex, '%e', A.nw);
fclose(findex);