Analysis of RP signals, profiles and LRC

The following script shows the basic data aquision from the Rake Probe

To get more detailed information about the Rake Probe visit its wikipage.

RP consists of one rake of Langmuir probes that can be operated in two basic modes.
The first one is floating potential mode ($V_\mathrm{float}$), where each pin is insulated
and we can measure an approximation of plasma potential.

The other mode is ion saturation current ($I_\mathrm{sat}$), where -100 V voltage is
applied to each pin. This voltage repels electrons so one can measure the plasma density.

Needed Python modules

import matplotlib.pyplot as plt
import numpy as np
import xarray as xr
try:
    import xrscipy.signal as dsp
except:
    print('xrscipy.signal module not found')
# import xarray_dsp.wavelets as wave

# for better readibility of spectrograms
from matplotlib.colors import LogNorm


# to keep the attributes
xr.set_options(keep_attrs=True)

# the data is downloaded in the directory/folder
# from where this script was started 
ds = np.DataSource('/tmp')

#inline plots
%matplotlib inline

Download (and save the data)

One might need to modify ulr_base variable in download_data function and the channel names in channel_list.

The correct url_base depends on the DAS used and the channel_list
represents the names of particular signals.
It is recommended to save the data to your computer if runnig the scipts locally,
as it saves time and makes your analysis independent on internet connection.

Eventually the data is converted into xarray.DataArray structure for convenience.

# test shot_no, the followin configuration should give the right data
# shot_no = 36443
# test channel_list
# channel_list = ['RP01', 'RP02','RP03','RP04','RP05','RP06', 'RP07',
#                 'RP08', 'RP09', 'RP10']
# test url_base
# "http://golem.fjfi.cvut.cz/shots/{}/Diagnostics/RakeProbe/".format(shot_no)

# a little helper function
def strip_off_data(base_url: str,data_name: str):
    """ This function is used to get data like parameters
    or comments from golem database"""
    url = base_url + data_name
    file = ds.open(url)
    parameters = np.loadtxt(file)
    return parameters

def get_basic_diagnostics(shot_no: int):
    url_base = f"http://golem.fjfi.cvut.cz/shots/{shot_no}/Diagnostics/BasicDiagnostics/Results/"
    
    basic_diagnostics = ["Bt", "Ip", "U_loop"]
    
    for diag in basic_diagnostics:...
        
    
def check_plasma_existence(shot_no: int, no_plasma_override: bool=False):
    """Checks, whether there was plasma, if there was not plasma
    raises RuntimeError. It is possible to override no_plasma and
    proceed in the code execution."""
    # base_url = "http://golem.fjfi.cvut.cz/shots/{}/Diagnostics/BasicDiagnostics/Results/"
    base_url = f"http://golem.fjfi.cvut.cz/shots/{shot_no}/Diagnostics/PlasmaDetection/Results/"    
    # check whether there was plasma
    
    b_plasma = strip_off_data(base_url,"b_plasma")
    if not b_plasma == 1.0:
        if not no_plasma_override:
            raise RuntimeError("There was no plasma.")
    
            


def download_data(shot_no: int, channel_list: list, no_plasma_override: bool=False):
    """Downloads the data and returns a tuple
    of time numpy.ndarray and data numpy.ndarray"""
    url_base = f"http://golem.fjfi.cvut.cz/shots/{shot_no}/Diagnostics/RakeProbe/"
    
    check_plasma_existence(shot_no, no_plasma_override)
    
    data_all = np.ones(shape=(40960,len(channel_list)))
    i = 0
    
    

    for channel in channel_list:
        # open the file
        file = ds.open(url_base+"{}.csv".format(channel))
        # transform to numpy format
        np_file = np.loadtxt(file, delimiter=',')
        # add to the complete array
        data_all[:,i] = np_file[:,1]
        # saves the time from the diag. output
        if i == 0:
            time = np_file[:,0] * 10e2 # to have time in ms
        i += 1
        
        # the 10e2 coefficient is given by a voltage divider
    return time, data_all * 1e2

# list of channel names in the database
channel_list = ['RP01', 'RP02','RP03','RP04','RP05','RP06', 'RP07',
                'RP08', 'RP09', 'RP10']
# shot is set to 0
# i.e. the last shot
# attention! when running the script locally
# leaving the 0 will not update the downloaded data to the current shot
# because it first looks wheather the data has already been downloaded

shot_no = 43652


time, data = download_data(shot_no, channel_list ,False)

# creation of xarray.DataArray
# it is with offset
data_all_off = xr.DataArray(data,
                            coords={'time':time,
                                    'signal':channel_list},
                            dims=['time','signal']
#                             ,attrs = {'shot_no'}
                            )

Download the time of the plasma start and the end

def get_plasma_start_end(shot_no: int, no_plasma_override: bool=False):
    """Downloads the plasma start and end if there was plasma."""
    # the following should not change in principle, but
    # if your are having problems with downloading the data,
    # check the url as adviced above
    base_url = f"http://golem.fjfi.cvut.cz/shots/{shot_no}/Diagnostics/PlasmaDetection/Results/"
    
    # checks whether there was plasma
    # and if yes, get the plasma start and end
    # this check is done here again because it is
    # meaningless to try to download plasma_start and plasma_end
    # if there was no plasma,
    check_plasma_existence(shot_no, no_plasma_override)
    
    t_start = strip_off_data(base_url,"t_plasma_start")
    t_end = strip_off_data(base_url,"t_plasma_end")
    return t_start, t_end

# add the plasma start and end as a attribute of 
# all_data xarray.DataArray
t_start, t_end = get_plasma_start_end(shot_no,False)


data_all_off.attrs = {"plasma_start":t_start, "plasma_end":t_end}

Susbtraction of the offset

# time window of period before the plasma start
# where we calculate the offset we subtract
time_off = slice(data_all_off.attrs['plasma_start']-2.5,data_all_off.attrs['plasma_start']-1)

# all the downloaded data without offset
data_all = data_all_off - data_all_off.sel(time=time_off).mean(dim='time')

Visual check of the signals

time_window = slice(data_all.attrs['plasma_start']-0.5,data_all.attrs['plasma_end']+1)
# the following list just to sort the pin coordinates such that
# it is shown as it is physically in the vessel
signal_order = ['RP01', 'RP02','RP03','RP04','RP05','RP06', 'RP07',
                'RP08', 'RP09', 'RP10']

# here is the part, where the plasma is is selected
sel_plasma = data_all.sel(time=time_window,signal=signal_order)


# plot and plot setup
all_plot = sel_plasma.plot.line(x="time", col='signal', col_wrap=2,
                                figsize=(12,12),
#                                 subplot_kws={'gridspec_kw:{'wspace':0.1}}
                               )

fontsize = 13
# axes labels
all_plot.set_xlabels('time [ms]',fontsize=fontsize)
all_plot.set_ylabels('$V_\mathrm{float}$ [V]',fontsize=fontsize)


# make the figure nicer
plt.tight_layout()

# save the figure with the name expected by DoubleRakeProbe.sh
plt.savefig('icon-fig.png')

Profile construction

To create a profile from the data, we need to choose the flat-top part of the discharge.
To do this, we find the maximum of the plasma current $I_{\mathrm{pl}_\mathrm{max}}$ and we find the largest time interval that $(t_\mathrm{start}, t_\mathrm{stop})$, that satisfies the following condition $$\left(\forall t \in \left(t_\mathrm{start}, t_\mathrm{stop}\right)\right)\left(I_\mathrm{pl}(t) \geq 0.95I_{\mathrm{pl}_\mathrm{max}}\right).$$ Then we perform ensemble average over $I_\mathrm{sat}(t)$, where $t \in \left(t_\mathrm{start}, t_\mathrm{stop}\right)$, for the time series we have this is equivalent to performing time average.

From the plotting point of view we need to know the position of the deepest Langmuir probe tip to have the correct values on the x axis.

# deepest_tip = strip_off_data(url_for_rake_probe_params, name_of_the_param)

# max_i_pl = 

Spectrograms of the data

we keep using the signal order and time window from the cell above → so we continue to work with the sel_plasma variable

spectrogram = dsp.spectrogram(sel_plasma, nperseg=512,dim='time')

spct_plot = spectrogram.plot(x='time',col='signal',col_wrap=2, norm=LogNorm(), figsize=(12,12), cbar_kwargs={'label':r'_[-]'} )

fontsize = 13
# axes labels
spct_plot.set_xlabels('time [ms]',fontsize=fontsize) spct_plot.set_ylabels('$f$ [kHz]',fontsize=fontsize)

plt.ylim(0,75)