def running_mean(l, N):
# Also works for the(strictly invalid) cases when N is even.
if (N//2)*2 == N:
N = N - 1
front = np.zeros(N//2)
back = np.zeros(N//2)
for i in range(1, (N//2)*2, 2):
front[i//2] = np.convolve(l[:i], np.ones((i,))/i, mode = 'valid')
for i in range(1, (N//2)*2, 2):
back[i//2] = np.convolve(l[-i:], np.ones((i,))/i, mode = 'valid')
return np.concatenate([front, np.convolve(l, np.ones((N,))/N, mode = 'valid'), back[::-1]])
def mkdir(dir):
try:os.makedirs(dir)
except OSError:pass
# integration of signal with the rid of offset
def integrate(array, dt = 1e-6, aver_line = 4500, coeff = 1):
aver = np.sum(array[0:aver_line])/aver_line
return np.cumsum(array - aver) * dt * coeff
def GetData(ShotNo,VacuumShot):
ShotNoStr=str(ShotNo)
os.system('wget http://golem.fjfi.cvut.cz/shots/' + ShotNoStr + '/basicdiagn/Btoroidal.npz -O '+ShotNoStr+'/Basics/Btoroidal.npz');
os.system('wget http://golem.fjfi.cvut.cz/shots/' + ShotNoStr + '/basicdiagn/Iplasma.npz -O '+ShotNoStr+'/Basics/Iplasma.npz');
os.system('wget http://golem.fjfi.cvut.cz/shots/' + ShotNoStr + '/basicdiagn/Uloop.npz -O '+ShotNoStr+'/Basics/Uloop.npz');
os.system('wget http://golem.fjfi.cvut.cz/shots/' + ShotNoStr + '/basicdiagn/graphpres.png -O '+ShotNoStr+'/Basics/graphres.png');
#Plasma parameters
PlasmaStart=int(float(np.loadtxt(urlopen(BaseURL+ShotNoStr+'/plasma_start')))*1e6) # in us
#PlasmaStart=10000
PlasmaEnd=int(float(np.loadtxt(urlopen(BaseURL+ShotNoStr+'/plasma_end')))*1e6) # in us
#PlasmaEnd=25000
np.savetxt(ShotNoStr+'/Basics/PlasmaStart',[PlasmaStart], delimiter=" ", fmt="%s")
np.savetxt(ShotNoStr+'/Basics/PlasmaEnd',[PlasmaEnd], delimiter=" ", fmt="%s")
os.system('wget '+BaseURL + ShotNoStr + '/' + CurrentTag + ' -O '+ShotNoStr+'/Basics/VAchar_current.txt');
os.system('wget '+BaseURL + ShotNoStr + '/' + VoltageTag + ' -O '+ShotNoStr+'/Basics/VAchar_voltage.txt');
os.system('wget '+BaseURL + ShotNoStr + '/' + PotentialTag + ' -O '+ShotNoStr+'/Basics/Potential.txt');
os.system('wget '+BaseURL + ShotNoStr + '/mirnov_5 -O '+ShotNoStr+'/Basics/mc5_p.txt');
os.system('wget '+BaseURL + str(VacuumShot) + '/mirnov_5 -O '+ShotNoStr+'/Basics/mc5_v.txt');
os.system('wget '+BaseURL + ShotNoStr + '/mirnov_13 -O '+ShotNoStr+'/Basics/mc13_p.txt');
os.system('wget '+BaseURL + str(VacuumShot) + '/mirnov_13 -O '+ShotNoStr+'/Basics/mc13_v.txt');
def DischargeBasics(ShotNo,VacuumShot,PlasmaStart, PlasmaEnd, Time, Uloop, Btoroidal, Iplasma, \
VAchar_voltage, VAchar_current,Potential, Potential_smoothed, dz):
f = plt.figure(figsize=(20.0, 5.0))
f,ax = plt.subplots(7,sharex=True);plt.subplots_adjust(hspace=0.001)
f.set_size_inches(FigSize)
ax[0].set_title('#' + str(ShotNo))
ax[0].plot(Time/1000,Uloop['data'][PlasmaStart:PlasmaEnd],label='Loop voltage');ax[0].set_ylabel('$U_l$ [V]');ax[0].set_ylim(0,);ax[0].legend(loc=0)
ax[1].plot(Time/1000,Btoroidal['data'][PlasmaStart:PlasmaEnd],label='Toroidal magnetic field');ax[1].set_ylabel('$B_t$ [T]');ax[1].set_ylim(0,);ax[1].legend(loc=0)
ax[2].plot(Time/1000,Iplasma['data'][PlasmaStart:PlasmaEnd]/1000,label='Plasma current');ax[2].set_ylabel('$I_p$ [kA]');ax[2].set_ylim(0,);ax[2].legend(loc=0)
ax[3].plot(Time/1000,VAchar_voltage[0:,1]*HWNapetovyDelic,label='VA char: U');ax[3].set_ylabel('$U$ [V]');ax[3].legend(loc=0)
ax[4].plot(Time/1000,VAchar_current[0:,1]/HWProudovyResistor*CurrentScale,label='VA char: I');ax[4].set_ylabel('$I$ [mA]');ax[4].legend(loc=0)
ax[4].plot(Time/1000,VAchar_current[0:,1]*0)
ax[5].plot(Time/1000,Potential[0:,1]*HWNapetovyDelic,label='Floating potential');ax[5].legend(loc=0)
ax[5].plot(Time/1000,Potential_smoothed*0)
ax[5].plot(Time/1000,Potential_smoothed,label='Floating potential smoothed');ax[5].set_ylabel('$<U_{fl}>$ [V]');ax[5].legend(loc=1)
ax[6].plot(Time/1000,dz[PlasmaStart:PlasmaEnd],label='Plasma position');ax[6].set_ylabel('PlPos [mm]');ax[6].legend(loc=0)
ax[6].set_xlabel('$t$ [ms]');#ax[6].set_ylim(-200,200)
ax[6].set_xlim(PlasmaStart/1000,PlasmaEnd/1000)
plt.savefig(str(ShotNo)+'/Basics/BasicDiagWithBPP.jpg', bbox_inches='tight')
if is_interactive():plt.show();
plt.close();
if is_interactive():plt.figure(figsize=(10, 8), dpi= 80, facecolor='w', edgecolor='k')
plt.plot(VAchar_voltage[:,1],VAchar_current[:,1],'.')
plt.xlabel('VA char - voltage U [V]');plt.ylabel('VA char - current U [V]')
plt.title('VA char raw XY data: current vs. voltage')
plt.savefig(str(ShotNo)+'/Basics/RawData.jpg', bbox_inches='tight')
if is_interactive():plt.show();
plt.close();
def FindingSweepIntervals(ShotNo,PlasmaStart, PlasmaEnd, Time, VAchar_voltage):
#ack http://stackoverflow.com/questions/4624970/finding-local-maxima-minima-with-numpy-in-a-1d-numpy-array
data=running_mean(VAchar_voltage[:,1]*HWNapetovyDelic,SmoothCoefficient)
Extrems_tmp = (np.diff(np.sign(np.diff(data))).nonzero()[0] + 1) # local min+max
# Sometimes double Extrems appear, removal:
#print Extrems
Extrems=[]
for i in range(len(Extrems_tmp)-1):
if abs(Extrems_tmp[i]-Extrems_tmp[i+1])>20 and abs(data[Extrems_tmp[i]])>60 :
Extrems.append(Extrems_tmp[i])
np.savetxt(str(ShotNo)+'/Data/Extrems',Extrems, fmt="%s")
# graphical output...
plt.figure(figsize=(FigSize), dpi= 80, facecolor='w', edgecolor='k')
plt.plot(Time/1000,data)
plt.plot(Time/1000,VAchar_voltage[:,1]*HWNapetovyDelic)
plt.plot(Time[Extrems]/1000, data[Extrems], "o", label="Extrems")
plt.ylabel('VAchar_voltage U [V]');plt.xlabel('time [us]')
plt.title('Maxima and minima identification')
plt.legend()
plt.savefig(str(ShotNo)+'/Basics/MaximMinims.jpg', bbox_inches='tight')
if is_interactive():plt.show();
plt.close();
return Extrems
def Fits(IndvSetup,PlasmaParameters,Extrems,Time,VAchar_voltage,VAchar_current,VAchar_voltage_smoothed,VAchar_current_smoothed, Potential_smoothed,Uloop, Btoroidal, Iplasma, dz, IndxDef):
ShotNo=IndvSetup[0]
ShotNoStr=str(ShotNo)
PlasmaStart=PlasmaParameters[0]
PlasmaEnd=PlasmaParameters[1]
ylim_min=np.min(VAchar_current_smoothed)
ylim_max=np.max(VAchar_current_smoothed)
xlim_min=np.min(VAchar_voltage_smoothed)
xlim_max=np.max(VAchar_voltage_smoothed)
print("Fits .. ...", end=":")
print(len(Extrems)-1)
AllVAchars=[]
def AzoozFn(Vbias, a1, a2, a3, a4):
return np.exp(a1 * np.tanh((Vbias+a2)/a3)) - a4
def ClassicalFn(Vbias, Isat, Vfl, Te,):
return -Isat*(1-np.exp((Vbias-Vfl)/Te))
def Classical4Fn(Vbias, Isat, Vfl, Te, C):
return -Isat*(1-C*(Vbias-Vfl))*(1-np.exp((Vbias-Vfl)/Te))
Results=[]
#Voltage = np.linspace(int(xlim_min),int(xlim_max),1000);dVoltage = Voltage[1]-Voltage[0]
#AzoozInitialGuess = [3.8,34,42,2.9]
AzoozInitialGuess = IndvSetup[3]
#AzoozInitialGuess = [1.5,8,24,1.15] # Absolutne nejslabsi misto
AFP0=0;AFP1=0;AFP2=0;AFP3=0
#ClassicalInitialGuess = [4.6,-30,21] # Absolutne nejslabsi misto
ClassicalInitialGuess = IndvSetup[4]
ClassicalFitParms=[0,0,0,0]
#ClassicalInitialGuess = [0.3,-20,20] # Absolutne nejslabsi misto
CFP0=0;CFP1=0;CFP2=0;CFP3=0;Pocet=0
if IndxDef[1]==0:
IndxDef[1]=len(Extrems)-1 # test purposes
for i in range(IndxDef[0],IndxDef[1]):
data=[VAchar_voltage_smoothed[Extrems[i]:Extrems[i+1]],\
VAchar_current_smoothed[Extrems[i]:Extrems[i+1]]]
#data=[VAchar_voltage[Extrems[i]:Extrems[i+1]],\
# VAchar_current[Extrems[i]:Extrems[i+1]]]
#Voltage = np.linspace(int(xlim_min),int(xlim_max),len(data[0]));dVoltage = Voltage[1]-Voltage[0]
Voltage = np.linspace(int(np.min(data[0])),int(np.max(data[0])),len(data[0]));dVoltage = Voltage[1]-Voltage[0]
#print(Voltage)
#Voltage = data[0];dVoltage = Voltage[1]-Voltage[0]
#print(Voltage)
def plotsetup():
plt.legend()
ax1.plot(data[0],data[1],',',label="rawdata")
ax1.plot(VAchar_voltage_smoothed[Extrems[i]:Extrems[i+1]],VAchar_current_smoothed[Extrems[i]:Extrems[i+1]],label="data averaged")
plt.grid(True)
ax1.set_xlim(xlim_min,xlim_max);
ax1.set_ylim(ylim_min,ylim_max*1.25); # Kvuli legende ..
ax1.set_title('VAchar #'+ShotNoStr+' \
t=<'+str(PlasmaStart+Extrems[i])+','+\
str(PlasmaStart+Extrems[i+1])+'> us')
ax1.set_xlabel('$U$ [V]')
ax1.set_ylabel('$I$ [mA]')
ax1.axhline(0);
ax2 = fig.add_axes([0.25, 0.45, 0.3, 0.18]);
ax2.set_xlim(PlasmaStart,PlasmaEnd);
ax2.set_ylim(0,6000);
ax2.set_yticks([j for j in xrange(0,6000,2000)])
ax2.set_yticklabels([j for j in xrange(0,6,2)])
ax2.set_xticks([j for j in xrange(PlasmaStart,PlasmaEnd,5000)])
ax2.set_xticklabels([j for j in xrange(PlasmaStart/1000,PlasmaEnd/1000,2)])
ax2.plot(Iplasma['data']);
ax2.axvline(PlasmaStart+(Extrems[i]+Extrems[i+1])/2);
ax2.set_ylabel('$I_p$ [kA]')
ax3 = fig.add_axes([0.25, 0.63, 0.3, 0.18]);
ax3.set_yticks([(j/10.0) for j in xrange(0,6,2)])
ax3.set_xticklabels([]);
ax3.set_ylim(0,0.55);
ax3.set_xlim(PlasmaStart,PlasmaEnd);
ax3.plot(Btoroidal['data']);
ax3.axvline(PlasmaStart+(Extrems[i]+Extrems[i+1])/2);
ax3.set_ylabel('$B_t$ [T]')
ax2.set_xlabel('$t$ [ms]')
#plt.show(); #just for tuning
fig, ax1 = plt.subplots();
plotsetup()
ax1.legend()
plt.savefig(ShotNoStr+'/IndivVAchars/VA'\
+str(i).zfill(2) +':'+str(PlasmaStart+Extrems[i])+'_'+str(PlasmaStart+Extrems[i+1])+'.jpg',\
bbox_inches='tight')
TimeLineDir=ShotNoStr+'/TimeLine/'+str(PlasmaStart+Extrems[i])+'_'+str(PlasmaStart+Extrems[i+1])
mkdir(TimeLineDir)
plt.savefig(TimeLineDir+'/VAcharSmoothed.jpg',bbox_inches='tight')
plt.close();
# Azooz fit
#if is_interactive():
print(i, end=",")
#else:
# print(i)
#print data[0][0]
AzoozFitParms, po_cov = spo.curve_fit(AzoozFn, data[0], data[1], AzoozInitialGuess)
#print(AzoozFitParms, po_cov)
AzoozFit=AzoozFn(Voltage,AzoozFitParms[0],AzoozFitParms[1],AzoozFitParms[2],AzoozFitParms[3])
#AzoozFit=AzoozFn(data[0],AzoozFitParms[0],AzoozFitParms[1],AzoozFitParms[2],AzoozFitParms[3])
dAzoozFit_dVoltage = np.gradient(AzoozFit, dVoltage)
InflPointLoc=argrelextrema(dAzoozFit_dVoltage, np.greater)[0]
#VflPointLoc=fsolve(AzoozFit, 0.0)
#VflPointLoc=fsolve(AzoozFn(Voltage,AzoozFitParms[0],AzoozFitParms[1],AzoozFitParms[2],AzoozFitParms[3]), 0.0)
#FitLimit=VflPointLoc
#print(len(InflPointLoc))
if len(InflPointLoc) != 1 \
or (ShotNo==23191 and i==7)\
or (ShotNo==23191 and i==17)\
or (ShotNo==23193 and i==11)\
or (ShotNo==23193 and i==13)\
or (ShotNo==23193 and i==19)\
or (ShotNo==23195 and i==7)\
or (ShotNo==23195 and i==9)\
or (ShotNo==23195 and i>22)\
or (ShotNo==23198 and i==7)\
or (ShotNo==23198 and i==8)\
or (ShotNo==23198 and i==12)\
or (ShotNo==23198 and i==13)\
or (ShotNo==23198 and i>22):\
# not reasonable data
print('(Azooz fit problem)',end="")
Results.append([\
i,\
(PlasmaStart+Extrems[i])/1000.0,\
(PlasmaStart+(Extrems[i]+Extrems[i+1])/2)/1000.0,\
(PlasmaStart+Extrems[i+1])/1000.0,\
nan,nan,nan,nan,nan,nan,nan,nan,nan,nan,nan,nan,nan,nan,nan,nan
])
np.savetxt(TimeLineDir+'/index.html',['\
<html><center>\
<h2><a href="../../index.html">Home</a></h2>\
<h1>\
<a href="../'+str(PlasmaStart+Extrems[i-1])+'_'+str(PlasmaStart+Extrems[i])+'/index.html"><=</a>\
'+str(i)+':Particular time <i>t=</i>'+str(Results[j][2])+' us\
<a href="../'+str(PlasmaStart+Extrems[i+1])+'_'+str(PlasmaStart+Extrems[i+2 if i<IndxDef[1]-1 else i])+'/index.html">=></a>\
</h1>\
<h2>Raw and smoothed data</h2>\
<a href="VAcharSmoothed.jpg"><img src="VAcharSmoothed.jpg" width="30%"/></a>\
<br/>!NO FITS!</center>'],fmt="%s")
else: # reasonable data
#nearest_idx = np.where(abs(data[0]-float(Voltage[InflPointLoc]))==abs(data[0]-float(Voltage[InflPointLoc])).min())[0][0] # Azooz limit
#nearest_idx = np.where(abs(data[0]-FitLimit)==abs(data[0]-FitLimit).min())[0][0] # HS limit
#nearest_idx = np.where(abs(AzoozFit-0.0)==abs(AzoozFit-0.0).min())[0][0] # Azooz V_fl limit
#print(nearest_idx)
#print(Voltage[nearest_idx])
errorflag=0
Shift=int(float(len(data[0]))/(np.max(data[0])-np.min(data[0]))*VoltageShift)
DataIndex=0
if data[0][0]<0: # Are we sweeping up or down??
#if 1==1:
#Voltage = np.linspace(int(np.max(data[0])),int(np.min(data[0])),len(data[0]));
nearest_idx = np.where(abs(AzoozFit-0.0)==abs(AzoozFit-0.0).min())[0][0] # Azooz V_fl limit
VoltageReduced=Voltage[0:nearest_idx+Shift]
DataIndex=int((nearest_idx+Shift)*1.0*len(data[0])/len(Voltage))
#print(nearest_idx,nearest_idx+Shift,DataIndex)
#print(VoltageReduced)
#VoltageReduced=VoltageDef[0:nearest_idx+Shift]
#print(VoltageReduced)
try:
#ClassicalFitParms, po_cov = spo.curve_fit(Classical4Fn,data[0][0:nearest_idx+Shift], data[1][0:nearest_idx+Shift], ClassicalInitialGuess)
ClassicalFitParms, po_cov = spo.curve_fit(Classical4Fn,data[0][0:DataIndex], data[1][0:DataIndex], ClassicalInitialGuess)
except RuntimeError:
errorflag=1
else:
#Voltage = np.linspace(int(np.min(data[0])),int(np.max(data[0])),len(data[0]));
#VoltageReduced=data[0][0:nearest_idx]
nearest_idx = np.where(abs(AzoozFit-0.0)==abs(AzoozFit-0.0).min())[0][0] # Azooz V_fl limit
#VoltageReduced=Voltage[0:len(data[0])-nearest_idx-Shift]
VoltageReduced=Voltage[0:nearest_idx+Shift]
DataIndex=int((len(data[0])-(nearest_idx+Shift))*1.0*len(data[0])/len(Voltage))
#print(VoltageReduced)
#VoltageReduced=VoltageDef[len(data[0])-nearest_idx-Shift:]
#print(VoltageReduced)
try:
#ClassicalFitParms, po_cov = spo.curve_fit(Classical4Fn,data[0][len(data[0])-nearest_idx-Shift:], data[1][len(data[0])-nearest_idx-Shift:], ClassicalInitialGuess)
ClassicalFitParms, po_cov = spo.curve_fit(Classical4Fn,data[0][DataIndex:], data[1][DataIndex:], ClassicalInitialGuess)
except RuntimeError:
errorflag=1
if errorflag==1:
print('(Classic fit problem)',end="")
Results.append([\
i,\
(PlasmaStart+Extrems[i])/1000.0,\
(PlasmaStart+(Extrems[i]+Extrems[i+1])/2)/1000.0,\
(PlasmaStart+Extrems[i+1])/1000.0,\
nan,nan,nan,nan,nan,nan,nan,nan,nan,nan,nan,nan,nan,nan,nan,nan
])
np.savetxt(TimeLineDir+'/index.html',['\
<html><center>\
<h2><a href="../../index.html">Home</a></h2>\
<h1>\
<a href="../'+str(PlasmaStart+Extrems[i-1])+'_'+str(PlasmaStart+Extrems[i])+'/index.html"><=</a>\
'+str(i)+':Particular time <i>t=</i>'+str(Results[j][2])+' us\
<a href="../'+str(PlasmaStart+Extrems[i+1])+'_'+str(PlasmaStart+Extrems[i+2 if i<IndxDef[1]-1 else i])+'/index.html">=></a>\
</h1>\
<h2>Raw and smoothed data</h2>\
<a href="VAcharSmoothed.jpg"><img src="VAcharSmoothed.jpg" width="30%"/></a>\
<br/>!NO FITS!</center>'],fmt="%s")
else:
ClassicalFit=Classical4Fn(VoltageReduced,ClassicalFitParms[0],ClassicalFitParms[1],ClassicalFitParms[2],ClassicalFitParms[3])
fig, ax1 = plt.subplots();
plotsetup()
ax1.plot(VoltageReduced, ClassicalFit, label="Classic fit")
ax1.plot(Voltage, AzoozFit,label="Azooz fit")
ax1.legend()
ax1.axvline(float(Voltage[InflPointLoc]))
ax1.axvline(float(Voltage[nearest_idx]))
ax1.axvline(float(Voltage[nearest_idx+Shift]))
plt.savefig(ShotNoStr+'/IndivVAcharsFits/VA'+str(i).zfill(2) +':'+str(PlasmaStart+Extrems[i])+'_'+str(PlasmaStart+Extrems[i+1])+'.jpg', bbox_inches='tight')
plt.savefig(ShotNoStr+'/TimeLine/'+str(PlasmaStart+Extrems[i])+'_'+str(PlasmaStart+Extrems[i+1])+'/VAcharFit.jpg', bbox_inches='tight')
np.savetxt(TimeLineDir+'/VoltageAzooz.dat',np.transpose([Voltage,AzoozFit]),fmt="%s")
np.savetxt(TimeLineDir+'/VoltageReducedClassic4.dat',np.transpose([VoltageReduced,ClassicalFit]),fmt="%s")
np.savetxt(TimeLineDir+'/SmoothedVachar.dat',np.transpose([VAchar_voltage_smoothed[Extrems[i]:Extrems[i+1]],VAchar_current_smoothed[Extrems[i]:Extrems[i+1]]]),fmt="%s")
np.savetxt(TimeLineDir+'/RawVachar.dat',np.transpose([VAchar_voltage[Extrems[i]:Extrems[i+1]],VAchar_current[Extrems[i]:Extrems[i+1]]]),fmt="%s")
np.savetxt(TimeLineDir+'/DataToBeClass4Fitted.dat',np.transpose([data[0][0:DataIndex], data[1][0:DataIndex]]),fmt="%s")
np.savetxt(TimeLineDir+'/nearest_idx.dat',[nearest_idx],fmt="%s")
np.savetxt(TimeLineDir+'/DataIndex.dat',[DataIndex],fmt="%s")
np.savetxt(TimeLineDir+'/Voltage@nearest_idx.dat',[Voltage[nearest_idx]],fmt="%s")
np.savetxt(TimeLineDir+'/Shift.dat',[Shift],fmt="%s")
np.savetxt(TimeLineDir+'/Voltage@nearest_idxPlusShift.dat',[Voltage[nearest_idx+Shift]],fmt="%s")
plt.close()
# plot Raw and smoothed data
plt.plot(VAchar_voltage[Extrems[i]:Extrems[i+1]],VAchar_current[Extrems[i]:Extrems[i+1]],'+', label="raw data")
plt.plot(VAchar_voltage_smoothed[Extrems[i]:Extrems[i+1]],VAchar_current_smoothed[Extrems[i]:Extrems[i+1]], label="smoothed")
plt.xlabel('$U$ [V]');plt.ylabel('$I$ [mA]');plt.legend()
plt.savefig(TimeLineDir+'/RawAndSmoothed.jpg', bbox_inches='tight')
plt.close()
# plot Raw and smoothed data with fits
plt.plot(VAchar_voltage[Extrems[i]:Extrems[i+1]],VAchar_current[Extrems[i]:Extrems[i+1]],',', label="raw data")
plt.plot(VAchar_voltage_smoothed[Extrems[i]:Extrems[i+1]],VAchar_current_smoothed[Extrems[i]:Extrems[i+1]], label="smoothed")
plt.plot(VoltageReduced, ClassicalFit, label="Classic fit")
plt.plot(Voltage, AzoozFit,label="Azooz fit")
plt.axvline(float(Voltage[nearest_idx]))
plt.axvline(float(Voltage[nearest_idx+Shift]))
plt.axhline(0)
#plt.plot(VAchar_voltage_smoothed[Extrems[i]:Extrems[i+1]],AzoozFn(data[0],AzoozFitParms[0],AzoozFitParms[1],AzoozFitParms[2],AzoozFitParms[3]), label="Azooz fit")
#plt.plot(VAchar_voltage_smoothed[Extrems[i]:Extrems[i+1]],Classical4Fn(data[0],ClassicalFitParms[0],ClassicalFitParms[1],ClassicalFitParms[2],ClassicalFitParms[3]), label="Classic4p fit")
plt.xlabel('$U$ [V]');plt.ylabel('$I$ [mA]');plt.legend()
plt.savefig(TimeLineDir+'/RawAndSmoothedWithFits.jpg', bbox_inches='tight')
plt.close()
# plot Raw and smoothed data with fits zoom
plt.plot(VAchar_voltage[Extrems[i]:Extrems[i+1]],VAchar_current[Extrems[i]:Extrems[i+1]],',', label="raw data")
plt.plot(VAchar_voltage_smoothed[Extrems[i]:Extrems[i+1]],VAchar_current_smoothed[Extrems[i]:Extrems[i+1]], label="smoothed")
plt.plot(VoltageReduced, ClassicalFit, label="Classic 4 fit")
plt.plot(Voltage, AzoozFit,label="Azooz fit")
plt.axvline(float(Voltage[nearest_idx]))
plt.axvline(float(Voltage[nearest_idx+Shift]))
plt.axhline(0)
plt.ylim(-10,10);
plt.xlabel('$U$ [V]');plt.ylabel('$I$ [mA]');plt.legend()
plt.savefig(TimeLineDir+'/RawAndSmoothedWithFitsZoom.jpg', bbox_inches='tight')
plt.close()
#plot Azooz stuff
plt.plot(VAchar_voltage_smoothed[Extrems[i]:Extrems[i+1]],AzoozFn(data[0],AzoozFitParms[0],AzoozFitParms[1],AzoozFitParms[2],AzoozFitParms[3]), label="Azooz fit")
plt.plot(Voltage,dAzoozFit_dVoltage*100,label="Azooz 1st derivative")
plt.xlabel('$U$ [V]');plt.ylabel('$I$ [mA]');plt.legend()
plt.savefig(TimeLineDir+'/AzoozStuff.jpg', bbox_inches='tight')
plt.close()
DischTime=(PlasmaStart+(Extrems[i]+Extrems[i+1])/2)/1000.0
V_BPP=sum(Potential_smoothed[(Extrems[i]):(Extrems[i+1])])/\
len(Potential_smoothed[(Extrems[i]):(Extrems[i+1])])
Results.append([\
i,\
(PlasmaStart+Extrems[i])/1000.0,\
DischTime,\
(PlasmaStart+Extrems[i+1])/1000.0,\
float(Btoroidal['data'][DischTime*1000:DischTime*1000+1]),\
AzoozFitParms[0],AzoozFitParms[1],AzoozFitParms[2],AzoozFitParms[3],\
#int(InflPointLoc),\
Voltage[nearest_idx],\
nearest_idx,\
float(Voltage[InflPointLoc]),\
ClassicalFitParms[0],ClassicalFitParms[1],ClassicalFitParms[2],ClassicalFitParms[3],\
(np.exp(AzoozFitParms[0])+AzoozFitParms[3])/(np.exp(-AzoozFitParms[0])+AzoozFitParms[3]),\
V_BPP,\
(V_BPP-ClassicalFitParms[1])/ClassicalFitParms[2],\
sum(dz[(PlasmaStart+Extrems[i]):(PlasmaStart+Extrems[i+1])])/\
len(dz[(PlasmaStart+Extrems[i]):(PlasmaStart+Extrems[i+1])])\
])
j=i-IndxDef[0] # Tuning purposes ..
k=0
np.savetxt(TimeLineDir+'/Results',Results[j],fmt="%s")
np.savetxt(TimeLineDir+'/index.html',['\
<html>\
<script type="text/x-mathjax-config">MathJax.Hub.Config({tex2jax: {inlineMath: [[\'$\',\'$\'], [\'\\(\',\'\\)\']]}});</script>\
<script src="https://cdn.mathjax.org/mathjax/latest/MathJax.js?config=TeX-AMS-MML_HTMLorMML"></script>\<body><center>\
<h2><a href="../../index.html">Home</a></h2>\
<h1>\
<a href="../'+str(PlasmaStart+Extrems[i-1])+'_'+str(PlasmaStart+Extrems[i])+'/index.html"><=</a>\
'+str(i)+':Particular time <i>t=</i>'+str(Results[j][2])+' us\
<a href="../'+str(PlasmaStart+Extrems[i+1])+'_'+str(PlasmaStart+Extrems[i+2 if i<IndxDef[1]-1 else i])+'/index.html">=></a>\
</h1>\
<h2>Raw and smoothed data with fits</h2>\
<a href="RawAndSmoothed.jpg"><img src="RawAndSmoothed.jpg" width="30%"/></a>\
<!--<h2>Raw and smoothed data with fits</h2>-->\
<a href="RawAndSmoothedWithFits.jpg"><img src="RawAndSmoothedWithFits.jpg" width="30%"/></a>\
<a href="RawAndSmoothedWithFitsZoom.jpg"><img src="RawAndSmoothedWithFitsZoom.jpg" width="30%"/></a><br/>\
<h2>Smoothed data with fits</h2>\
<a href="VAcharFit.jpg"><img src="VAcharFit.jpg" width="30%"/></a><br/>\
<h2>Results</h2>\
<b>i</b>:'+'{:2d}'.format(Results[j][0]) +\
', <b>Time:</b>From:'+'{:4.2f}'.format(Results[j][1]) +\
',via:'+'{:4.2f}'.format(Results[j][2]) +\
',to:'+'{:4.2f}'.format(Results[j][3]) +\
',$B_t$:'+'{:4.2f}'.format(Results[j][4]) +\
'</br><b>Azooz fit:</b> $a_1$='+'{:+5.3e}'.format(Results[j][5]) +\
',$a_2$='+'{:4.2f}'.format(Results[j][6]) +\
',$a_3$='+'{:+5.3e}'.format(Results[j][7]) +\
',$a_4$='+'{:4.2f}'.format(Results[j][8]) +\
'</br><b>Inflex point:</b>'+'{:+4.0f}'.format(Results[j][9]) +\
',Infl'+'{:+4.0f}'.format(Results[j][10]) +\
',Infl'+'{:+3.2f}'.format(Results[j][11]) +\
'</br><b>Class fit:</b> $I_{sat}$='+'{:4.2f}'.format(Results[j][12]) +\
',$\\varphi_{pl}$='+'{:4.2f}'.format(Results[j][13]) +\
',$T_e$='+'{:4.2f}'.format(Results[j][14]) +\
',$C$='+'{:4.2f}'.format(Results[j][15]) +\
',$Ratio_{sat}$='+'{:4.2f}'.format(Results[j][16]) +\
',$\\varphi_{fl}$='+'{:4.2f}'.format(Results[j][17]) +\
',$\\alpha$='+'{:4.2f}'.format(Results[j][18]) +\
',Plasma position='+'{:4.2f}'.format(Results[j][19]) +\
'<h2>Azooz stuff</h2>\
<br/><a href="AzoozStuff.jpg"><img src="AzoozStuff.jpg" width="30%"/></a><br/>\
<h2>Data</h2>\
<a href="Results">Results</a><br/>\
<a href="../'+str(PlasmaStart+Extrems[i-1])+'_'+str(PlasmaStart+Extrems[i])+'/">Characteristics</a><br/>\
</center></body></html>'],fmt="%s")
AFP0=AFP0+AzoozFitParms[0];AFP1=AFP1+AzoozFitParms[1];AFP2=AFP2+AzoozFitParms[2];AFP3=AFP3+AzoozFitParms[3]
CFP0=CFP0+ClassicalFitParms[0];CFP1=CFP1+ClassicalFitParms[1];CFP2=CFP2+ClassicalFitParms[2];CFP3=CFP3+ClassicalFitParms[3];Pocet=Pocet+1
#print(end=",")
#np.savetxt(TimeLineDir+'/chars.dat',\
# [[VAchar_voltage[Extrems[i]:Extrems[i+1]][j],\
# VAchar_current[Extrems[i]:Extrems[i+1]][j],\
# VAchar_voltage_smoothed[Extrems[i]:Extrems[i+1]][j],\
# VAchar_current_smoothed[Extrems[i]:Extrems[i+1]][j],\
# AzoozFn(data[0][j],AzoozFitParms[0],AzoozFitParms[1],AzoozFitParms[2],AzoozFitParms[3]), \
# dAzoozFit_dVoltage[int(1000*(j*1.0)/len(data[0]))],\
# Classical4Fn(data[0][j],ClassicalFitParms[0],ClassicalFitParms[1],ClassicalFitParms[2],ClassicalFitParms[3])] \
# for j in xrange(len(data[0]))],delimiter=" ",fmt="%.3f")
np.savetxt(ShotNoStr+'/Data/Results',Results,fmt="%+8.2f")
#print(IndvSetup[3][0])
np.savetxt(ShotNoStr+'/Data/AzoozFitCoefficients',['Fit parameters generation *estimate*-><final>: No:{:d}. <br/>'.format(Pocet)\
+'Azooz fit: a1:{:4.2f} -> {:4.2f};'.format(IndvSetup[3][0],AFP0/Pocet)\
+'a2:{:4.2f} -> {:4.2f};'.format(IndvSetup[3][1],AFP1/Pocet)\
+'a3:{:4.2f} -> {:4.2f};'.format(IndvSetup[3][2],AFP2/Pocet)\
+'a4:{:4.2f} -> {:4.2f};<br/>'.format(IndvSetup[3][3],AFP3/Pocet)\
+'Classic fit: Isat:{:4.2f} -> {:4.2f};'.format(IndvSetup[4][0],CFP0/Pocet)\
+'varphi_pl:{:4.2f} -> {:4.2f};'.format(IndvSetup[4][1],CFP1/Pocet)\
+'Te:{:4.2f} -> {:4.2f};'.format(IndvSetup[4][2],CFP2/Pocet)\
+'C:{:4.3f} -> {:4.3f};'.format(IndvSetup[4][3],CFP3/Pocet)\
],fmt="%s")
return Results
def htmlgeneration(ShotNo,UBt,VacuumShot,PlasmaStart, PlasmaEnd, Time, Extrems, IndxDef, Results):
print("html generation ...")
ShotNoStr=str(ShotNo)
os.system('rm '+str(ShotNo)+'/index.html');
fileid = open(str(ShotNo)+'/index.html','a+')
fileid.write('<html><head><title>VA char @ GOLEM</title>\
<style>\
div.grafy30 {display:inline-block;width: 30%;}\
div.grafy20 {display:inline-block;width: 20%;}\
.grafy20 img {width: 100%;}\
.grafy30 img {width: 100%;}\
</style>\
<meta http-equiv="Content-Type" content="text/html; charset=utf-8">\
<style></style>\
<script type="text/x-mathjax-config">MathJax.Hub.Config({tex2jax: {inlineMath: [[\'$\',\'$\'], [\'\\(\',\'\\)\']]}});</script>\
<script src="https://cdn.mathjax.org/mathjax/latest/MathJax.js?config=TeX-AMS-MML_HTMLorMML"></script>\
</head><body><center>')
fileid.write('<h2><a href="../index.html">Home</a></h2>')
fileid.write('<h1>BPP experiments @ the tokamak GOLEM</h1>')
fileid.write('<h3>Essential parameters:</h3><table><tr><td><ul>\
<li>Current resistor $R=$'+str(HWProudovyResistor)+' $\Omega$</li>\
<li>Voltage divider: 1:'+str(HWNapetovyDelic)+'</li>\
<li>The shot: <a href="http://golem.fjfi.cvut.cz/shots/'+ShotNoStr+'">'+ShotNoStr+'</a></li>\
<li>The DAS used: <a href="http://golem.fjfi.cvut.cz/shots/'+ShotNoStr+'/DAS/1011Papouch_St.ON/">Papouch_St</a></li>\
<li>$U_{B_t}$='+str(UBt)+'</li>\
<li>Vacuum shot /plasma position issues/:<a href="http://golem.fjfi.cvut.cz/shots/'+str(VacuumShot)+'">'+str(VacuumShot)+'</a></li>\
<li>Voltage shift='+str(VoltageShift)+'</li>\
</ul></td></tr></table>')
fileid.write('<div class="grafy30"><h2>Experimental setup</h2>')
fileid.write('<h3>Ball pen probe @ North-East port</h3>')
fileid.write('<center>\
<a href="http://golem.fjfi.cvut.cz/wikiraw/Experiments/EdgePlasmaPhysics/ParticleFlux/BallPenProbe/Experiments/BottomPosition/setup/vI/ExpSetup-BPP.png">\
<img src="../setup/ExpSetup-BPP.png" width="100%"></a><br>\
<a href="/wikiraw/Experiments/EdgePlasmaPhysics/ParticleFlux/BallPenProbe/Experiments/BottomPosition/setup/">setup</a></div>')
fileid.write('<div class="grafy30"><h2>The discharge #'+ShotNoStr+'</h2>')
fileid.write('<h3>Basic diagnostics</h3>\
<a href="http://golem.fjfi.cvut.cz/shots/'+ShotNoStr+'/"><img src="Basics/graphres.png" width="100%"></a><br/>\
Plasma start:<a href="ShotNo/PlasmaStart">'+str(PlasmaStart)+'</a> us\
Plasma end:<a href="ShotNo/PlasmaEnd">'+str(PlasmaEnd)+'</a> us</div>\
<h3>The BPP data vs diagnostics</h3>\
<a href="Basics/BasicDiagWithBPP.jpg"><img src="Basics/BasicDiagWithBPP.jpg" width="30%"></a><br/>\
<a href="ShotNo/">Reference shot data</a><br/>\
<div class="grafy30"><h3>The BPP raw data</h3>\
<a href="ShotNo/RawData.jpg"><img src="Basics/RawData.jpg"></a><br/></div>')
fileid.write('\
<div class="grafy30"><h3>Maxima and minima localization in bpp_voltage</h3>\
<a href="ShotNo/MaximMinims.jpg"><img src="Basics/MaximMinims.jpg"></a><br/>\
<a href="Data/Extrems">Extrems</a><br/>\
<a href="Data/minima">minima</a><br/>\
<a href="Data/maxima">maxima</a></div>\
')
fileid.write('<h3>Final movie</h3>')
fileid.write('<img src="finalmovie.gif"><br/><a href="IndivVAchars/">Individual figures</a></br>\
.dat format: VAchar_voltage,VAchar_current,VAchar_voltage_smoothed,VAchar_current_smoothed,Azooz fit,dAzoozFit/dVoltage,ClassicFit<br/><div>')
StartIndx=0
if IndxDef[1]==0:
IndxDef[1]=len(Extrems)-1 # test purposes
for i in range(IndxDef[0]-IndxDef[0],IndxDef[1]-IndxDef[0]):
if os.path.isfile(str(ShotNo)+'/IndivVAcharsFits/VA'+str(i).zfill(2) +':'+str(PlasmaStart+Extrems[i])+'_'+str(PlasmaStart+Extrems[i+1])+'.jpg'): # Azooz fit problem
fileid.write('<div class="grafy20">\
<center><a href="TimeLine/'+str(PlasmaStart+Extrems[i])+'_'+str(PlasmaStart+Extrems[i+1])+'/index.html">\
<img src="IndivVAcharsFits/VA'+str(i).zfill(2) +':'+str(PlasmaStart+Extrems[i])+'_'+str(PlasmaStart+Extrems[i+1])+'.jpg" width="100%"></a>\
</br>\
<b>i</b>:'+'{:2d}'.format(Results[i][0]) +\
', <b>Time:</b>From:'+'{:4.2f}'.format(Results[i][1]) +\
',via:'+'{:4.2f}'.format(Results[i][2]) +\
',to:'+'{:4.2f}'.format(Results[i][3]) +\
',$B_t$'+'{:4.2f}'.format(Results[i][4]) +\
'</br><b>Azooz fit:</b> $a_1$='+'{:+5.3e}'.format(Results[i][5]) +\
',$a_2$='+'{:4.2f}'.format(Results[i][6]) +\
',$a_3$='+'{:+5.3e}'.format(Results[i][7]) +\
',$a_4$='+'{:4.2f}'.format(Results[i][8]) +\
'</br><b>Class fit:</b> $I_{sat}$='+'{:4.2f}'.format(Results[i][12]) +\
',$\\varphi_{pl}$='+'{:4.2f}'.format(Results[i][13]) +\
',$T_e$='+'{:4.2f}'.format(Results[i][14]) +\
',$Ratio_{sat}$='+'{:4.2f}'.format(Results[i][15]) +\
',$\\varphi_{fl}$='+'{:4.2f}'.format(Results[i][16]) +\
',$\\alpha$='+'{:4.2f}'.format(Results[i][17]) +\
'<br><a href=TimeLine/'+str(PlasmaStart+Extrems[i])+'_'+str(PlasmaStart+Extrems[i+1])+'/chars.dat>data</a></center></div>')
else:
fileid.write('<div class="grafy20">\
<center><a href="TimeLine/'+str(PlasmaStart+Extrems[i])+'_'+str(PlasmaStart+Extrems[i+1])+'/index.html">\
<img src="IndivVAchars/VA'+str(i).zfill(2) +':'+str(PlasmaStart+Extrems[i])+'_'+str(PlasmaStart+Extrems[i+1])+'.jpg" width="100%"></a>\
</br>\
<b>i</b>:'+'{:2d}'.format(Results[i][0]) +\
', <b>Time:</b>From:'+'{:4.2f}'.format(Results[i][1]) +\
',via:'+'{:4.2f}'.format(Results[i][2]) +\
',to:'+'{:4.2f}'.format(Results[i][3]) +\
'<br/><b>!!Azooz fit problem! -> No Fits !</b><br/><br/><br/><br/>' +\
'<br><a href=Data/VA'+str(i).zfill(2) +':'+str(PlasmaStart+Extrems[i])+'_'+str(PlasmaStart+Extrems[i+1])+'.dat>data</a></center></div>')
fileid.write('</div>\
<img src="GlobalResults.jpg"></a><br/>')
fileid.write('<table border="1">\
<tr><th>$i$</th><th colspan=3>Time</th><th>$B_t$</th>\
<th colspan=4>Azooz fit parameters<br/>$I(a,V_{bias})=\exp[a_1 \\tanh\{(V_{bias}+a_2)/a_3\}]+a_4$</th>\
<th colspan=3>InflexPoint</th>\
<th colspan=4>Classical fit parameters<br/>$I(V_{bias},I_{sat},V_{fl},T_e,C)=-I_{sat}*(1-C*(V_{bias}-V_{fl}))*(1-\exp((V_{bias}-V_{fl})/T_e)$</th>\
<th>$Ratio_{sat}$</th><th>$\\varphi_{fl}$</th><th>$\\alpha$</th><th>$Pl_{pos}$ [mm]</th></tr>\
<tr><th></th><th>Start</th><th>Middle</th><th>End</th>\
<th></th>\
<th>$a_1$</th><th>$a_2$</th><th>$a_3$</th><th>$a_4$</th>\
<th>Data row</th><th>Data row</th><th>Voltage</th>\
<th>$I_{sat}$</th><th>$\\varphi_{pl}$</th><th>$T_e$</th><th>C</th><th></th><th></th><th></th><th></th></tr>')
for i in range(IndxDef[0]-IndxDef[0],IndxDef[1]-IndxDef[0]):
fileid.write('<tr>\
<td>'+'{:2.0f}'.format(Results[i][0]) +'</td>\
<td>'+'{:4.2f}'.format(Results[i][1]) +'</td>\
<td><a href=TimeLine/'+str(PlasmaStart+Extrems[i])+'_'+str(PlasmaStart+Extrems[i+1])+'/index.html>'+'{:4.2f}'.format(Results[i][2]) +'</a></td>\
<td>'+'{:4.2f}'.format(Results[i][3]) +'</td>\
<td>'+'{:4.2f}'.format(Results[i][4]) +'</td>\
<td>'+'{:4.2f}'.format(Results[i][5]) +'</td>\
<td>'+'{:4.2f}'.format(Results[i][6]) +'</td>\
<td>'+'{:4.2f}'.format(Results[i][7]) +'</td>\
<td>'+'{:4.2f}'.format(Results[i][8]) +'</td>\
<td>'+'{:3.0f}'.format(Results[i][9]) +'</td>\
<td>'+'{:3.0f}'.format(Results[i][10]) +'</td>\
<td>'+'{:4.2f}'.format(Results[i][11]) +'</td>\
<td>'+'{:4.2f}'.format(Results[i][12]) +'</td>\
<td>'+'{:4.2f}'.format(Results[i][13]) +'</td>\
<td>'+'{:4.2f}'.format(Results[i][14]) +'</td>\
<td>'+'{:4.3f}'.format(Results[i][15]) +'</td>\
<td>'+'{:4.2f}'.format(Results[i][16]) +'</td>\
<td>'+'{:4.2f}'.format(Results[i][17]) +'</td>\
<td>'+'{:4.2f}'.format(Results[i][18]) +'</td>\
<td>'+'{:4.2f}'.format(Results[i][19]) +'</td>\
</tr>')
fileid.write('</table><br/><a href="Data/Results">data</a><br/>')
fileid.write(np.loadtxt(ShotNoStr+'/Data/AzoozFitCoefficients',dtype='str'))
fileid.write('</div><h2>Experimental Photo</h2>')
fileid.write('<center><img src="http://golem.fjfi.cvut.cz/wikiraw/Experiments/EdgePlasmaPhysics/ParticleFlux/BallPenProbe/Experiments/BottomPosition/setup/vI/2017-01-19_23-11-12_w.jpg"></center>\n\n')
fileid.write('<h2>Resources</h2>')
fileid.write('<ul>\
<li><a href="VAcharSweeped.ipynb">Jupyter notebook</a><li>\
<li><a href="VAcharSweeped.py">Jupyter notebook python export</a><li>\
<li><a href="VAcharSweeped.html">Jupyter notebook html export</a><li>\
<li><a href="VAcharSweeped.pdf">Jupyter notebook pdf export</a><li>\
</ul>')
fileid.write('<h2>References</h2>')
fileid.write('<ul>\
<li><a href="https://www.mathworks.com/matlabcentral/fileexchange/19312-langmuir-probe-data-analysis-code">Aasim Azooz: Langmuir probe data analysis code. Mathworks.</a></li>\
</ul>')
fileid.write('</body></html>')
fileid.close()
def GlobalGraphGen(ShotNo,VacuumShot,PlasmaStart, PlasmaEnd, Time, Extrems, Uloop, Btoroidal, Iplasma,\
VAchar_current, Potential, Potential_smoothed, PlasmaPosition, dz, Results):
print("Global graph generation .. ")
ShotNoStr=str(ShotNo)
Borders=0,len(Extrems) # Initial and final problems out (in future option)
ResultsArr=np.array(Results)
f = plt.figure(figsize=(20.0, 10.0))
f,ax = plt.subplots(17,sharex=True);plt.subplots_adjust(hspace=0.001)
f.set_size_inches(10,20)
ax[0].set_title('#' + ShotNoStr)
ax[0].plot(Time/1000,Uloop['data'][PlasmaStart:PlasmaEnd]);ax[0].set_ylabel('$U_l$ [V]')
ax[1].plot(Time/1000,Btoroidal['data'][PlasmaStart:PlasmaEnd]);ax[1].set_ylabel('$B_t$ [T]')
ax[1].plot(ResultsArr[:,2],ResultsArr[:,4],'+')
ax[2].plot(Time/1000,Iplasma['data'][PlasmaStart:PlasmaEnd]/1000);ax[2].set_ylabel('$I_p$ [kA]')
ax[3].plot(Time/1000,VAchar_current[0:,1]/HWProudovyResistor*CurrentScale,label='VA char: I');ax[3].set_ylabel('$I$ [mA]');ax[3].legend(loc=0)
ax[3].plot(Time/1000,VAchar_current[0:,1]*0)
ax[4].plot(Time/1000,dz[PlasmaStart:PlasmaEnd]);ax[4].set_ylabel('$<Pl_{pos}>$ [mm]')
ax[4].plot(Time/1000,dz[PlasmaStart:PlasmaEnd]*0)
ax[4].plot(ResultsArr[:,2],ResultsArr[:,19],'+')
ax[5].plot(ResultsArr[:,2],ResultsArr[:,5],'+');ax[5].set_ylabel('$a_1$')
ax[6].plot(ResultsArr[:,2],ResultsArr[:,6],'+');ax[6].set_ylabel('$a_2$')
ax[7].plot(ResultsArr[:,2],ResultsArr[:,7],'+');ax[7].set_ylabel('$a_3$')
ax[8].plot(ResultsArr[:,2],ResultsArr[:,8],'+');ax[8].set_ylabel('$a_4$')
ax[9].plot(ResultsArr[:,2],ResultsArr[:,9],'+');ax[9].set_ylabel('InflPoint [V]')
ax[10].plot(ResultsArr[Borders[0]:Borders[1],2],ResultsArr[Borders[0]:Borders[1],12],'+');ax[10].set_ylabel('$I_{sat}$')
ax[11].plot(ResultsArr[Borders[0]:Borders[1],2],ResultsArr[Borders[0]:Borders[1],13],'+');ax[11].set_ylabel('$\\varphi_{pl}$ [V]')
ax[12].plot(ResultsArr[Borders[0]:Borders[1],2],ResultsArr[Borders[0]:Borders[1],14],'+');ax[12].set_ylabel('$T_e$ [eV]')
ax[13].plot(ResultsArr[Borders[0]:Borders[1],2],ResultsArr[Borders[0]:Borders[1],15],'+');ax[13].set_ylabel('$C$ []')
ax[14].plot(ResultsArr[Borders[0]:Borders[1],2],ResultsArr[Borders[0]:Borders[1],16],'+');ax[14].set_ylabel('$Ratio_{sat}$ []')
ax[15].plot(ResultsArr[Borders[0]:Borders[1],2],ResultsArr[Borders[0]:Borders[1],18],'+');ax[15].set_ylabel('$\\alpha$ []')
ax[16].plot(Time/1000,Potential_smoothed);ax[16].set_ylabel('$\\varphi_{fl}$ [V]')
ax[16].plot(Time/1000,Potential_smoothed*0)
ax[16].plot(ResultsArr[:,2],ResultsArr[:,17],'+')
ax[16].set_xlim(PlasmaStart/1000*8/10,PlasmaEnd/1000*11/10);ax[16].set_xlabel('$t$ [ms]')
plt.savefig(ShotNoStr+'/GlobalResults.jpg', bbox_inches='tight')
if is_interactive():plt.show();
plt.close();
def Initials(ShotNo):
ShotNoStr=str(ShotNo)
mkdir(ShotNoStr+'/IndivVAchars');mkdir(ShotNoStr+'/IndivVAcharsFits');
mkdir(ShotNoStr+'/Basics');mkdir(ShotNoStr+'/Data');mkdir(ShotNoStr+'/TimeLine');
def PlasmaParametersDef(ShotNo):
return [\
int(np.loadtxt(str(ShotNo)+'/Basics/PlasmaStart')),\
int(np.loadtxt(str(ShotNo)+'/Basics/PlasmaEnd'))\
]
def PlasmaPositionFn(ShotNo):
ShotNoStr=str(ShotNo)
mc5_p=np.loadtxt(ShotNoStr+'/Basics/mc5_p.txt')
mc5_v=np.loadtxt(ShotNoStr+'/Basics/mc5_v.txt')
mc13_p=np.loadtxt(ShotNoStr+'/Basics/mc13_p.txt')
mc13_v=np.loadtxt(ShotNoStr+'/Basics/mc13_v.txt')
#upper mirnov coil
mc5_p = integrate(mc5_p[0:,1], coeff=1/(3705 * 1e-6))
mc5_v = integrate(mc5_v[0:,1], coeff=1/(3705 * 1e-6))
mc5 = mc5_p - mc5_v # pure signal from plasma
#bottom mirnov coil
mc13_p = integrate(mc13_p[0:,1], coeff=1/(3705 * 1e-6))
mc13_v = integrate(mc13_v[0:,1], coeff=1/(3705 * 1e-6))
mc13 = mc13_p - mc13_v # pure signal from plasma
dz= 93 * (mc5 - mc13)/(mc5 + mc13) # calculation of displacement
for i in range(len(dz)):
if dz[i]>50:dz[i]=50
if dz[i]<-50:dz[i]=-50
return dz
def IndividualShot(IndvSetup,IndxDef):
ShotNo=IndvSetup[0];print( 'Doing', ShotNo)
VacuumShot=IndvSetup[1]
ShotNoStr=str(ShotNo)
Initials(ShotNo)
PlasmaParameters=PlasmaParametersDef(ShotNo)
PlasmaStart=PlasmaParameters[0]
PlasmaEnd=PlasmaParameters[1]
Time = np.linspace(PlasmaStart,PlasmaEnd,(PlasmaEnd-PlasmaStart))
VAchar_voltage=np.loadtxt(ShotNoStr+'/Basics/VAchar_voltage.txt')[PlasmaStart:PlasmaEnd]
VAchar_current=np.loadtxt(ShotNoStr+'/Basics/VAchar_current.txt')[PlasmaStart:PlasmaEnd]
VAchar_voltage_smoothed=running_mean(VAchar_voltage[:,1]*HWNapetovyDelic,SmoothCoefficient)
VAchar_current_smoothed=running_mean(VAchar_current[:,1]/HWProudovyResistor*CurrentScale,SmoothCoefficient)
Btoroidal=np.load(ShotNoStr+'/Basics/Btoroidal.npz');
Iplasma=np.load(ShotNoStr+'/Basics/Iplasma.npz');
Uloop=np.load(ShotNoStr+'/Basics/Uloop.npz');
Potential=np.loadtxt(ShotNoStr+'/Basics/Potential.txt')[PlasmaStart:PlasmaEnd]
Potential_smoothed=running_mean(Potential[:,1]*HWNapetovyDelic,SmoothCoefficient)
PlasmaPosition=PlasmaPositionFn(ShotNo)
DischargeBasics(ShotNo,VacuumShot,PlasmaStart, PlasmaEnd, Time, Uloop, Btoroidal, Iplasma,\
VAchar_voltage, VAchar_current,Potential, Potential_smoothed, PlasmaPosition)
Extrems=FindingSweepIntervals(ShotNo,PlasmaStart, PlasmaEnd, Time, VAchar_voltage)
Results=Fits(IndvSetup,PlasmaParameters,Extrems,Time,VAchar_voltage[:,1]*HWNapetovyDelic,VAchar_current[:,1]/HWProudovyResistor*CurrentScale,VAchar_voltage_smoothed,VAchar_current_smoothed,\
Potential_smoothed, Uloop, Btoroidal, Iplasma, PlasmaPosition, IndxDef)
htmlgeneration(ShotNo,IndvSetup[2],VacuumShot,PlasmaStart, PlasmaEnd, Time, Extrems, IndxDef, Results)
GlobalGraphGen(ShotNo,VacuumShot,PlasmaStart, PlasmaEnd, Time, Extrems, Uloop, Btoroidal, Iplasma,\
VAchar_current,Potential, Potential_smoothed, PlasmaPosition, PlasmaPosition, Results)
def GlobalFigs():
mkdir('GlobalFigs')
Time = np.linspace(0,40000,40000)
for graph in GraphSetup:
FigDef()
for shot in Setup:
data=np.load(str(shot[0])+'/Basics/'+graph[0]+'.npz')
plt.plot(Time/1000,data['data'],label='#'+str(shot[0])+' ($U_{B_t}$='+str(shot[2])+')');
plt.legend(loc=0)
plt.title('Reproducibility for: '+graph[0])
plt.xlabel('$t$ [ms]')
plt.ylabel(graph[1])
plt.savefig('GlobalFigs/Reproducibility_'+graph[0]+'.jpg', bbox_inches='tight')
if is_interactive():plt.show();
plt.close();
def Dependences():
# Dependences with respect to time (Results[0:,2])
for result in ResultsDef:
FigDef()
for shot in Setup:
Results=np.loadtxt(str(shot[0])+'/Data/Results')
plt.plot(Results[0:,2],Results[0:,result[0]],shot[6],label=str(shot[0])+' ($U_{B_t}$='+str(shot[2])+')')
plt.title('Crosscheck for: '+result[2])
plt.ylabel(result[2])
plt.xlabel('$t$ [ms]')
plt.legend(loc=0)
plt.savefig('GlobalFigs/TimeCrossCheck_'+str(result[1])+'.jpg', bbox_inches='tight')
if is_interactive():plt.show();
plt.close();
os.system('rm GlobalFigs/TimeCrossCheck_'+str(result[1])+'.dat')
with open('GlobalFigs/TimeCrossCheck_'+str(result[1])+'.dat', "a") as myfile:
for i in xrange(len(Results[0:,3])):
k=1;datarow=str(Results[i][3])
#print(lst)
for shot in Setup:
datarow=datarow+result[4].format(Results[0:,result[0]][i])
k=k+1
myfile.write(datarow+'\n')
# and Ufl x Plpot speciality:
FigDef()
for shot in Setup:
Results=np.loadtxt(str(shot[0])+'/Data/Results')
plt.plot(Results[0:,2],Results[0:,13],shot[6],label=str(shot[0])+' $\\varphi_{fl}$')
plt.plot(Results[0:,2],Results[0:,17],shot[6],label=str(shot[0])+' $\\varphi_{pl}$')
plt.title('Crosscheck for: $\\varphi_{fl}$ x $\\varphi_{pl}$')
plt.ylabel('$V_{fl}$ [V]')
plt.xlabel('$t$ [ms]')
plt.legend(loc=0)
plt.savefig('GlobalFigs/TimeCrossCheck_FlowPot_PlPot.jpg', bbox_inches='tight')
if is_interactive():plt.show();
plt.close();
# Dependences with respect to Bt (Results[0:,4])
for result in ResultsDef:
FigDef()
for shot in Setup:
Results=np.loadtxt(str(shot[0])+'/Data/Results')
plt.plot(Results[0:,4],Results[0:,result[0]],shot[6],label=str(shot[0])+' ($U_{B_t}$='+str(shot[2])+')')
plt.title('Crosscheck for: '+result[2])
plt.ylabel(result[2])
plt.xlabel('$B_t$ [T]')
plt.legend(loc=0)
plt.savefig('GlobalFigs/BtCrossCheck_'+str(result[1])+'.jpg', bbox_inches='tight')
if is_interactive():plt.show();
plt.close();
os.system('rm GlobalFigs/BtCrossCheck_'+str(result[1])+'.dat')
with open('GlobalFigs/BtCrossCheck_'+str(result[1])+'.dat', "a") as myfile:
for i in xrange(len(Results[0:,4])):
k=1;datarow=str(Results[i][4])
for shot in Setup:
datarow=datarow+result[4].format(Results[0:,result[0]][i])
k=k+1
myfile.write(datarow+'\n')
# and Ufl x Plpot speciality:
FigDef()
for shot in Setup:
Results=np.loadtxt(str(shot[0])+'/Data/Results')
plt.plot(Results[0:,4],Results[0:,13],shot[6],label=str(shot[0])+' $\\varphi_{fl}$')
plt.plot(Results[0:,4],Results[0:,17],shot[6],label=str(shot[0])+' $\\varphi_{pl}$')
plt.title('Crosscheck for: $\\varphi_{fl}$ x $\\varphi_{pl}$')
plt.ylabel('$V_{fl}$ [V]')
plt.xlabel('$B_t$ [T]')
plt.legend(loc=0)
plt.savefig('GlobalFigs/BtCrossCheck_FlowPot_PlPot.jpg', bbox_inches='tight')
if is_interactive():plt.show();
plt.close();
def HtmlGeneration():
print("html generation ...")
os.system('rm index.html');
fileid = open('index.html','a+')
fileid.write('<html><head><title>VA char @ GOLEM</title>\
<style>\
div.grafy30 {display:inline-block;width: 30%;}\
div.grafy20 {display:inline-block;width: 20%;}\
.grafy20 img {width: 100%;}\
.grafy30 img {width: 100%;}\
</style>\
<meta http-equiv="Content-Type" content="text/html; charset=utf-8">\
<style></style>\
<script type="text/x-mathjax-config">MathJax.Hub.Config({tex2jax: {inlineMath: [[\'$\',\'$\'], [\'\\(\',\'\\)\']]}});</script>\
<script src="https://cdn.mathjax.org/mathjax/latest/MathJax.js?config=TeX-AMS-MML_HTMLorMML"></script>\
</head><body><center>')
fileid.write('<h2>BPP experiments @ the tokamak GOLEM</h2>')
fileid.write('<h3>Essential parameters:</h3><table><tr><td><ul>\
<li>Date: 090217</li>\
<li>BPP used for Plasma Float potential $\\varphi_{fl}$, BPP_LP used in sweped mode</li>\
<li>Probe position:'+Probe_position+'</li>\
<li>The VA characteristic resistor:'+str(HWProudovyResistor)+'$\Omega$</li>\
<li>Sweeping frequency:'+str(HFSweepFrequency)+'</li>\
</ul></td></tr></table>')
fileid.write('<div class="grafy30"><h2>Experimental setup</h2>')
fileid.write('<h3>Ball pen probe @ North-East port</h3>')
fileid.write('<center>\
<a href="http://golem.fjfi.cvut.cz/wikiraw/Experiments/EdgePlasmaPhysics/ParticleFlux/BallPenProbe/Experiments/BottomPosition/setup/vI/ExpSetup-BPP.png">\
<img src="setup/ExpSetup-BPP.png" width="100%"></a><br>\
<a href="'+Base+'/setup/">setup</a></div>\
<h2>Discharges involved</h2>\
<table border="1"><tr>\
<th>Discharge #</th>\
<th>Vacuum #</th>\
<th>$U_{B_t}$ [V]</th>\
<th>Analysis</th>\
<th>Comment</th></tr>')
for shot in Setup:
fileid.write('<tr>\
<td><a href="http://golem.fjfi.cvut.cz/shots/'+str(shot[0])+'">#'+str(shot[0])+'</td>\
<td><a href="http://golem.fjfi.cvut.cz/shots/'+str(shot[1])+'">#'+str(shot[1])+'</td>\
<td>'+str(shot[2])+'</td>\
<td><a href="'+str(shot[0])+'/index.html">link</td>\
<td>'+str(shot[5])+'</td>\
</tr>')
fileid.write('</table><h2>Reproducibility</h2>')
for graph in GraphSetup:
fileid.write('<div class="grafy30"><h2>'+graph[0]+':</h2>\
<a href="GlobalFigs/Reproducibility_'+graph[0]+'.jpg">\
<img src="GlobalFigs/Reproducibility_'+graph[0]+'.jpg" width="100%"></a></div>')
fileid.write('<h2>Fits definition</h2>\
Azooz fit definition:$I(a,V_{bias})=\exp[a_1 \\tanh\{(V_{bias}+a_2)/a_3\}]+a_4$<br/>\
Classical fit definition:$I(V_{bias},I_{sat},V_{fl},T_e,C)=-I_{sat}*(1-C*(V_{bias}-V_{fl}))*(1-\exp((V_{bias}-V_{fl})/T_e)$')
fileid.write('<h2>Time cross compares</h2>')
for result in ResultsDef:
fileid.write('<div class="grafy30"><h2>'+str(result[1])+'</h2>\
<a href="GlobalFigs/TimeCrossCheck_'+str(result[1])+'.jpg">\
<img src="GlobalFigs/TimeCrossCheck_'+str(result[1])+'.jpg" width="100%"></a></br>\
<a href="GlobalFigs/TimeCrossCheck_'+str(result[1])+'.dat">data</a></div>')
fileid.write('<div class="grafy30"><h2> $\\varphi_{fl}$ x $\\varphi_{pl}$</h2>\
<a href="GlobalFigs/TimeCrossCheck_FlowPot_PlPot.jpg">\
<img src="GlobalFigs/TimeCrossCheck_FlowPot_PlPot.jpg" width="100%"></a></br></div>')
fileid.write('<h2>$B_t$ cross compares</h2>')
for result in ResultsDef:
fileid.write('<div class="grafy30"><h2>'+str(result[1])+'</h2>\
<a href="GlobalFigs/BtCrossCheck_'+str(result[1])+'.jpg">\
<img src="GlobalFigs/BtCrossCheck_'+str(result[1])+'.jpg" width="100%"></a></br>\
<a href="GlobalFigs/BtCrossCheck_'+str(result[1])+'.dat">data</a></div>')
fileid.write('<div class="grafy30"><h2> $\\varphi_{fl}$ x $\\varphi_{pl}$</h2>\
<a href="GlobalFigs/BtCrossCheck_FlowPot_PlPot.jpg">\
<img src="GlobalFigs/BtCrossCheck_FlowPot_PlPot.jpg" width="100%"></a></br></div>')
fileid.write('<h2>References</h2>')
fileid.write('<ul>\
<li><a href="https://www.mathworks.com/matlabcentral/fileexchange/19312-langmuir-probe-data-analysis-code">Aasim Azooz: Langmuir probe data analysis code. Mathworks.</a></li>\
</ul>')
fileid.write('</body></html>')
fileid.close()
def MovieGeneration(ShotNo):
os.system('convert -delay 50 '+str(ShotNo)+'/IndivVAchars/* '+str(ShotNo)+'/finalmovie.gif')