Skip to content

Commit 436cbe2

Browse files
ASvyatkovskiyASvyatkovskiy
committed
Initial version of MDSplus scripts from Kyle
1 parent f7a48c9 commit 436cbe2

File tree

2 files changed

+310
-0
lines changed

2 files changed

+310
-0
lines changed

‎data/get_data.py‎

Lines changed: 86 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,86 @@
1+
import MDSplus
2+
import numpy
3+
import time
4+
import sys
5+
6+
class gadata:
7+
"""GA Data Obj"""
8+
def __init__(self,signal,shot,tree=None,connection=None,nomds=False):
9+
10+
# Save object values
11+
self.signal = signal
12+
self.shot = shot
13+
self.zdata = -1
14+
self.xdata = -1
15+
self.ydata = -1
16+
self.zunits = ''
17+
self.xunits = ''
18+
self.yunits = ''
19+
self.rank = -1
20+
self.connection = connection
21+
22+
## Retrieve Data
23+
t0 = time.time()
24+
found = 0
25+
26+
# Create the MDSplus connection (thin) if not passed in
27+
if self.connection is None:
28+
self.connection = MDSplus.Connection('atlas.gat.com')
29+
30+
# Retrieve data from MDSplus (thin)
31+
if nomds == False:
32+
try:
33+
if tree != None:
34+
tag = self.signal
35+
fstree = tree
36+
else:
37+
tag = self.connection.get('findsig("'+self.signal+'",_fstree)').value
38+
fstree = self.connection.get('_fstree').value
39+
self.connection.openTree(fstree,shot)
40+
self.zdata = self.connection.get('_s = '+tag).data()
41+
self.zunits = self.connection.get('units_of(_s)').data()
42+
self.rank = numpy.ndim(self.zdata)
43+
self.xdata = self.connection.get('dim_of(_s)').data()
44+
self.xunits = self.connection.get('units_of(dim_of(_s))').data()
45+
if self.xunits == '' or self.xunits == ' ':
46+
self.xunits = self.connection.get('units(dim_of(_s))').data()
47+
if self.rank > 1:
48+
self.ydata = self.connection.get('dim_of(_s,1)').data()
49+
self.yunits = self.connection.get('units_of(dim_of(_s,1))').data()
50+
if self.yunits == '' or self.yunits == ' ':
51+
self.yunits = self.connection.get('units(dim_of(_s,1))').data()
52+
found = 1
53+
54+
# MDSplus seems to return 2-D arrays transposed. Change them back.
55+
if numpy.ndim(self.zdata) == 2: self.zdata = numpy.transpose(self.zdata)
56+
if numpy.ndim(self.ydata) == 2: self.ydata = numpy.transpose(self.ydata)
57+
if numpy.ndim(self.xdata) == 2: self.xdata = numpy.transpose(self.xdata)
58+
59+
except Exception,e:
60+
# print ' Signal not in MDSplus: %s' % (signal,)
61+
pass
62+
63+
64+
# Retrieve data from PTDATA
65+
if found == 0:
66+
self.zdata = self.connection.get('_s = ptdata2("'+signal+'",'+str(shot)+')')
67+
if len(self.zdata) != 1:
68+
self.xdata = self.connection.get('dim_of(_s)')
69+
self.rank = 1
70+
found = 1
71+
72+
# Retrieve data from Pseudo-pointname
73+
if found == 0:
74+
# print ' Signal not in PTDATA: %s' % (signal,)
75+
self.zdata = self.connection.get('_s = pseudo("'+signal+'",'+str(shot)+')')
76+
if len(self.zdata) != 1:
77+
self.xdata = self.connection.get('dim_of(_s)')
78+
self.rank = 1
79+
found = 1
80+
81+
if found == 0:
82+
# print ' Signal not in pseudo-pointname: %s' % (signal,)
83+
print " No such signal: %s" % (signal,)
84+
return
85+
86+
print ' GADATA Retrieval Time : ',time.time() - t0

‎data/get_mdsplus_data.py‎

Lines changed: 224 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,224 @@
1+
'''
2+
http://www.mdsplus.org/index.php?title=Documentation:Tutorial:RemoteAccess&open=76203664636339686324830207&page=Documentation%2FThe+MDSplus+tutorial%2FRemote+data+access+in+MDSplus
3+
http://piscope.psfc.mit.edu/index.php/MDSplus_%26_python#Simple_example_of_reading_MDSplus_data
4+
http://www.mdsplus.org/documentation/beginners/expressions.shtml
5+
http://www.mdsplus.org/index.php?title=Documentation:Tutorial:MdsObjects&open=76203664636339686324830207&page=Documentation%2FThe+MDSplus+tutorial%2FThe+Object+Oriented+interface+of+MDSPlus
6+
'''
7+
from __future__ import print_function
8+
from MDSplus import *
9+
from data_processing import ShotList
10+
from pylab import *
11+
import numpy as np
12+
import sys
13+
import multiprocessing as mp
14+
from functools import partial
15+
import Queue
16+
import os
17+
import errno
18+
19+
import gadata
20+
21+
#print("Importing numpy version"+np.__version__)
22+
23+
prepath = '/home/felkerk/mdsplus_data/'
24+
shot_numbers_path = 'shot_lists/'
25+
save_path = 'signal_data/'
26+
machine = 'd3d'
27+
28+
29+
if machine == 'nstx':
30+
shot_numbers_files = ['disrupt_nstx.txt']
31+
server_path = "skylark.pppl.gov:8501::"
32+
signal_paths = ['engineering/ip1/',
33+
'operations/rwmef_plas_n1_amp_br/',
34+
'efit02/li/',
35+
'activespec/ts_ld/',
36+
'passivespec/bolom_totpwr/',
37+
'nbi/nb_p_inj/',
38+
'efit02/wpdot/']
39+
40+
elif machine == 'd3d':
41+
# shot_numbers_files = ['shotlist_JaysonBarr_clear.txt']
42+
# shot_numbers_files = ['shotlist_JaysonBarr_disrupt.txt']
43+
shot_numbers_files = ['d3d_single_clear.txt']# ,'d3d_clear.txt', 'd3d_disrupt.txt']
44+
server_path = 'atlas.gat.com'
45+
# signal_paths = ['PINJ','IP','Q95','DENSITY','WMHD'] #,'PRAD'] #PRAD returns a 2D xdata
46+
# Recommended signals from DIII-D
47+
signal_paths = ['efsli','ipsip','efsbetan','efswmhd','nssampn1l','nssfrqn1l',
48+
'nssampn2l','nssfrqn2l',
49+
'dusbradial','dssdenest',r'\bol_l15_p',r'\bol_l03_p','bmspinj','bmstinj','pcechpwrf']
50+
elif machine == 'jet':
51+
shot_numbers_files = ['CWall_clear.txt','CFC_unint.txt','BeWall_clear.txt','ILW_unint.txt']
52+
server_path = 'mdsplus.jet.efda.org'
53+
54+
#plasma current, locked mode, output power
55+
signal_paths = ['jpf/da/c2-ipla',
56+
'jpf/da/c2-loca',
57+
'jpf/db/b5r-ptot>out']
58+
59+
60+
61+
#internal inductance, time derivative of stored energy, input power, total diamagnetic energy
62+
signal_paths += ['jpf/gs/bl-li<s',
63+
'jpf/gs/bl-fdwdt<s',
64+
'jpf/gs/bl-ptot<s',
65+
'jpf/gs/bl-wmhd<s']
66+
67+
68+
69+
70+
#density signals
71+
#4 vertical channels and 4 horizontal channels
72+
signal_paths += ['jpf/df/g1r-lid:{:03d}'.format(i) for i in range(1,9)]
73+
74+
75+
76+
#radiation signals
77+
#vertical signals, don't use signal 16 and 23
78+
signal_paths += ['jpf/db/b5vr-pbol:{:03d}'.format(i) for i in range(1,28) if (i != 16 and i != 23)]
79+
signal_paths += ['jpf/db/b5hr-pbol:{:03d}'.format(i) for i in range(1,24)]
80+
81+
82+
#ece temperature profiles
83+
#temperature of channel i vs time
84+
signal_paths += ['ppf/kk3/te{:02d}'.format(i) for i in range(1,97)]
85+
86+
#radial position of channel i vs time
87+
#signal_paths += ['ppf/kk3/ra{:02d}'.format(i) for i in range(1,97)]
88+
89+
#radial position of channel i mapped onto midplane vs time
90+
signal_paths += ['ppf/kk3/rc{:02d}'.format(i) for i in range(1,97)]
91+
92+
93+
94+
95+
96+
97+
98+
else:
99+
print('unkown machine. exiting')
100+
exit(1)
101+
102+
103+
104+
105+
106+
107+
def create_missing_value_filler():
108+
time = np.linspace(0,100,1000)
109+
vals = np.zeros_like(time)
110+
return time,vals
111+
112+
def mkdirdepth(filename):
113+
folder=os.path.dirname(filename)
114+
if not os.path.exists(folder):
115+
os.makedirs(folder)
116+
117+
118+
def get_tree_and_tag(path):
119+
spl = path.split('/')
120+
tree = spl[0]
121+
tag = '\\' + spl[1]
122+
return tree,tag
123+
124+
125+
def format_save_path(prepath,signal_path,shot_num):
126+
return prepath + signal_path + '/{}.txt'.format(shot_num)
127+
128+
129+
def save_shot(shot_num_queue,c,signal_paths,save_prepath,machine):
130+
missing_values = 0
131+
# if machine == 'd3d':
132+
# reload(gadata) #reloads Gadata object with connection
133+
while True:
134+
try:
135+
shot_num = shot_num_queue.get(False)
136+
except Queue.Empty:
137+
break
138+
for signal_path in signal_paths:
139+
save_path_full = format_save_path(save_prepath,signal_path,shot_num)
140+
if os.path.isfile(save_path_full):
141+
print('-',end='')
142+
else:
143+
try:
144+
if machine == 'nstx':
145+
tree,tag = get_tree_and_tag(signal_path)
146+
c.openTree(tree,shot_num)
147+
data = c.get(tag).data()
148+
time = c.get('dim_of('+tag+')').data()
149+
elif machine == 'jet':
150+
try:
151+
data = c.get('_sig=jet("{}/",{})'.format(signal_path,shot_num)).data()
152+
time = c.get('_sig=dim_of(jet("{}/",{}))'.format(signal_path,shot_num)).data()
153+
except:
154+
missing_values += 1
155+
print('Signal {}, shot {} missing. Filling with zeros'.format(signal_path,shot_num))
156+
time,data = create_missing_value_filler()
157+
elif machine == 'd3d':
158+
try:
159+
ga1 = gadata.gadata('{}'.format(signal_path),shot_num,tree='d3d',connection=c)
160+
# ga1 = gadata.gadata('\\{}'.format(signal_path),shot_num,tree='d3d',connection=c)
161+
data = ga1.zdata
162+
time = ga1.xdata
163+
except:
164+
missing_values += 1
165+
print('Signal {}, shot {} missing. Filling with zeros'.format(signal_path,shot_num))
166+
time,data = create_missing_value_filler()
167+
168+
data_two_column = vstack((time,data)).transpose()
169+
try: #can lead to race condition
170+
mkdirdepth(save_path_full)
171+
except OSError, e:
172+
if e.errno == errno.EEXIST:
173+
# File exists, and it's a directory, another process beat us to creating this dir, that's OK.
174+
pass
175+
else:
176+
# Our target dir exists as a file, or different error, reraise the error!
177+
raise
178+
savetxt(save_path_full,data_two_column,fmt = '%f %f')
179+
print('.',end='')
180+
except:
181+
print('Could not save shot {}, signal {}'.format(shot_num,signal_path))
182+
print('Warning: Incomplete!!!')
183+
raise
184+
sys.stdout.flush()
185+
print('saved shot {}'.format(shot_num))
186+
print('Finished with {} missing values total'.format(missing_values))
187+
188+
189+
190+
191+
192+
save_prepath = prepath+save_path + '/' + machine + '/'
193+
194+
shot_numbers,_ = ShotList.get_multiple_shots_and_disruption_times(prepath + shot_numbers_path,shot_numbers_files)
195+
196+
197+
fn = partial(save_shot,signal_paths=signal_paths,save_prepath=save_prepath,machine=machine)
198+
num_cores = min(mp.cpu_count(),8) #can only handle 8 connections at once :(
199+
queue = mp.Queue()
200+
for shot_num in shot_numbers:
201+
queue.put(shot_num)
202+
connections = [Connection(server_path) for _ in range(num_cores)]
203+
processes = [mp.Process(target=fn,args=(queue,connections[i])) for i in range(num_cores)]
204+
205+
print('running in parallel on {} processes'.format(num_cores))
206+
start_time = time.time()
207+
208+
for p in processes:
209+
p.start()
210+
for p in processes:
211+
p.join()
212+
213+
print('Finished downloading {} shots in {} seconds'.format(len(shot_numbers),time.time()-start_time))
214+
215+
# c = Connection(server_path)
216+
217+
# pool = mp.Pool()
218+
# for shot_num in shot_numbers:
219+
# # save_shot(shot_num,signal_paths,save_prepath,machine,c)
220+
# for (i,_) in enumerate(pool.imap_unordered(fn,shot_numbers)):
221+
# print('{}/{}'.format(i,len(shot_numbers)))
222+
223+
# pool.close()
224+
# pool.join()

0 commit comments

Comments
 (0)