#!/usr/bin/env python 
# -*- coding: utf-8 -*-



"""
Program for controling of the Ocean optics spectrometers and extractin the most important properties from spectra


"""



from numpy import *
import matplotlib
matplotlib.use('agg') 

from matplotlib.pyplot import *

import os
from time import mktime, strptime,sleep,asctime, localtime,struct_time
import re
from scipy.optimize import leastsq   
from scipy.interpolate import interpolate
from scipy import stats
import  scipy.stats 
from scipy.special import erf
import ConfigParser
import zipfile
import socket
import  shutil 
import locale
import subprocess

from matplotlib.ticker import NullFormatter,ScalarFormatter

set_printoptions(precision=4,linewidth=400)

rcParams['backend'] = 'Agg'




def MAD(x):
    """
    Median absolute deviation
    """
    return median(abs(x-median(x)))*1.4826
    
def MovingAverage(data, subset_size):
    """
    Simple moving average
    """
    if subset_size < 1 and len(data) < subset_size:
	raise ValueError('subset_size must be > 1 and less then len(data)')
    
    mvavg = zeros(size(data))
    for i in range(subset_size):
	mvavg[subset_size/2:-subset_size/2]+= data[i: (+i-subset_size)]
    
    mvavg/= float(subset_size)
    mvavg[:subset_size/2] = mean(data[:subset_size/2+1])
    mvavg[-subset_size/2:] = mean(data[-subset_size/2-1:])
    
    return mvavg
    
    
def UniversalLineShape( p,x):
    """
    The most general line shape. xc is position, A amplitude, w width, a skewness, n voitig profile
    Line profile is integrated over single pixels 
    """
    x = linspace(x[0],x[-1],size(x)*10) #increase resolution
    xc =p[2]
    A = p[0]
    w = abs(p[1])
    a = p[3]
    n = p[4]
    # "skew voitig profile"
    #http://en.wikipedia.org/wiki/Skew_normal_distribution
    y  = abs(A)*(1/sqrt(2*pi*w**2)*exp(-((x-xc)/w)**2/2)*(1+erf(a*(x-xc)/sqrt(2)/w))*(1-n)+ n*2/pi*w/(4*(x-xc)**2+w**2))
    y = reshape(y, (size(x)/10,10) ) 
    y = sum(y,1)/10 #integrate over pixel
    return y

def lineCharacteristics(p):
    """
    Calculate real properties of the line shape from function UniversalLineShape. 
    pc[0] - surface under line, pc[1] - FWHM, pc[2] - center of mass, pc[3] - skewness, pc[4] - voitig profile
    """
    pc = zeros(5)
    a = p[3]
    n = p[4]
    p[1] = abs(p[1])
    d = a/sqrt(1+a**2)
    pc[0] = p[0]
    pc[2] = p[2]+(1-n)*p[1]*d*sqrt(2/pi)
    f_L = p[1]/2*n
    f_G = (1-n)*p[1]*sqrt(1-2*d**2/pi)*2*sqrt(2*log(2))
    
    pc[1] = 0.5*f_L+sqrt(0.22*f_L**2+f_G**2)#http://en.wikipedia.org/wiki/Voigt_profile
    
    pc[3] = (1-n)*(4-pi)/2*(d*sqrt(2/pi))**3/(1-2*d**2/pi)**(1.5)
    pc[4] = n
    return pc
    


    
class Tokamak:
    "Tokamak properties and methods container class " 

    
    def __init__(self, name):
	self.name = name.upper()
	self.loadConfig()
	
	self.shotNumber = -1
	self.uploadServerIp =  '0.0.0.0'
	
	self.elements = list()

	for name in self.ElementsList:
	    self.elements.append([loadtxt('SpectralLines/'+name+'.txt'),name])
	
	
    
    def loadConfig(self):
	"""
	Load tokamak properties from config file
	"""
	if not os.path.isfile('tokamak_'+self.name+'.cfg'):
	    raise Exception('Tokamak config file does not exists.')
	
	config = ConfigParser.RawConfigParser()
	config.read('tokamak_'+self.name+'.cfg')
	name = config.get('DEFAULT', 'Name')
	if name != self.name:
	    raise Exception('WrongConfig')
	self.dataAcqusitionTime = config.getfloat('Basic parameters', 'Plasma length')
	self.trigger_shift = config.getfloat('Basic parameters', 'trigger shift')

	self.database = config.get('Database', 'remote database folder')
	self.trigger_port = config.getint('Database', 'trigger port')
	
	self.login = config.get('Database', 'database login')
	elements_str = config.get('Elements', 'ions')
	self.ElementsList = elements_str.replace("', '", ' ').strip("'[]").split()

	


    def saveConfig(self):
	"""
	Save tokamak properties to config file
	"""
	config = ConfigParser.RawConfigParser()
	config.set('DEFAULT', 'Name', self.name)
	
	config.add_section('Basic parameters')
	config.set('Basic parameters', 'Plasma length', self.dataAcqusitionTime)
	config.set('Basic parameters', 'trigger shift', self.trigger_shift)

	config.add_section('Database')
	config.set('Database', 'remote database folder', self.database)
	config.set('Database', 'trigger port', self.trigger_port)

	config.set('Database', 'database login', self.login)
	
	config.add_section('Elements')
	config.set('Elements', 'ions', self.ElementsList)

	
	with open('tokamak_'+self.name +'.cfg', 'wb') as configfile:
	    config.write(configfile)
	    
	
	    
	 
 ##GOLEM
 
    def storeData(self, data_folder):  
	"""
	Upload measured data and graphs to database
	"""
	#TODO odstranit!!!
	#self.uploadServerIp = '192.168.2.125'
	#self.shotNumber = 6879
	
	
	if self.uploadServerIp == '0.0.0.0' or self.shotNumber == -1: 
	    raise Exception('can not access to database')
	if os.path.isfile(data_folder+'spectra.txt_LinesEvolution_.png'):
	    os.system('convert -resize 150 %sspectra.txt_LinesEvolution_.png %sicon.png' %((data_folder,)*2))
	    
	
	while(True):
	    try:		
		if os.system(('scp -r %s %s:%s' % (data_folder,self.login+'@'+self.uploadServerIp,self.database ))%(self.shotNumber)):
			raise Exception('can not copy data')

		if os.path.isfile(data_folder+'icon.png'):
		    if os.system(('scp -r %s %s:%s' % (data_folder+'icon.png',self.login+'@'+self.uploadServerIp,self.database+'icon.png'))%(self.shotNumber)):
			raise Exception('can not copy icon')    
		break
		
	    except Exception as inst:
		print 'problem with database access:', inst 
		sleep(5)
		
		
    def actualShotNumber(self):  
    	"""
	Return number of the actual shot
	"""
	return self.shotNumber
	
    def waitOnSoftwareTrigger(self): 
    	"""
	Waits on the signal from the main server, that database is ready. Also shot number of actual shot  is reclieved. 
	"""
	def netcat(hostname, port):
	    s = socket.socket(socket.AF_INET, socket.SOCK_STREAM)
	    s.bind((hostname, port))
	    s.listen(1)
	    conn, addr = s.accept()
	    data = conn.recv(1024)
	    conn.close()
	    return addr,data
      
	HOST = ''                 # Symbolic name meaning all available interfaces
	PORT =  self.trigger_port            # Arbitrary non-privileged port
	
	#self.uploadServerIp = '192.168.2.125'   #TODO vrátit!!
	#self.shotNumber = 7341
	#return    
	[addr,data] = netcat(HOST, PORT)
	self.shotNumber = int(data)   # je to integer
	self.uploadServerIp = addr[0]
	
	
	
    def waitOnHardwareTrigger(self, filename):
	"""
	Detect that spectrometer has reclieved hardware trigger
	"""
	#filename = './loopnumber'
	#detect the hardware trigger only by change of the file with the data
	if not os.path.isfile(filename):
	    creat_time = 0
	else:
	    creat_time = os.stat(filename).st_mtime

	change_time = creat_time
	while(creat_time == change_time):
	    if not os.path.isfile(filename):
		change_time = 0
	    else:
		change_time = os.stat(filename).st_mtime
	    sleep(1)
	    
	
#class Tokamak:
    #def __init__(self, name):
	#self.name = name.upper()
	#self.loadConfig()
	
	#self.shotNumber = -1
	#self.uploadServerIp =  '0.0.0.0'
	
    
    #def loadConfig(self):
	
	#if not os.path.isfile('tokamak_'+self.name+'.cfg'):
	    #raise Exception('Tokamak config file does not exists.')
	
	#config = ConfigParser.RawConfigParser()
	#config.read('tokamak_'+self.name+'.cfg')
	#name = config.get('DEFAULT', 'Name')
	#if name != self.name:
	    #raise Exception('WrongConfig')
	#self.dataAcqusitionTime = config.getfloat('Basic parameters', 'Plasma length')
	#self.database = config.get('Database', 'remote database folder')
	#self.trigger_port = config.getint('Database', 'trigger port')
	#self.login = config.get('Database', 'database login')
	

	
	
	
	#TODO sezna očekávaných prvků! v golemovi není B, i v pořadí v jakém se nají hledat

	
	    
	 
 #COMPASS
 
    #def storeData(self, data_folder):  
    	"""
	Upload measured data and graphs to database
	"""
	#newpath = '../'+ str(self.shotNumber)+'/'
	#if os.path.exists(newpath):
	    #print 'file '+newpath+'exists, will be deleted'
	    #shutil.rmtree(newpath)
	    
	#shutil.copytree(data_folder,newpath)
	
		
		
    #def actualShotNumber(self):
	"""
	Return number of the actual shot
	"""
	#self.shotNumber = -1
	
	#if os.path.exists('index.php'):
	    #os.remove('index.php')
	
	#if os.system("wget --post-data='"+self.login.decode('base64')+ "' " +self.database) or not os.path.exists('index.php'):
	    #print 'Problems with logbook access'
	    #return self.shotNumber
	    
	    
	#f = open('index.php', 'r')
	#search_str = "name=shot_no class=input_text_s value="
	
	#for i,line in enumerate(f):
	    #if line.find(search_str)   != -1:
		#ind = line.find(search_str)  
		#self.shotNumber = int(line[ind+len(search_str): ind+4+len(search_str)])
		#break

	#if self.shotNumber == -1:
	    #raise  Exception('Problems with logbook access')	
	
	#return self.shotNumber
	
	
	
    #def waitOnSoftwareTrigger(self):
	"""
	Waits on the signal from the main server, that database is ready. Also shot number of actual shot  is reclieved. 
	"""
	#pass
	
	
	
    #def waitOnHardwareTrigger(self, filename):
	"""
	Detect that spectrometer has reclieved hardware trigger
	"""
	##detect the hardware trigger only by change of the file with the data
	#if not os.path.isfile(filename):
	    #creat_time = 0
	#else:
	    #creat_time = os.stat(filename).st_mtime

	#change_time = creat_time
	#while(creat_time == change_time):
	    #if not os.path.isfile(filename):
		#change_time = 0
	    #else:
		#change_time = os.stat(filename).st_mtime
	    #sleep(1)
	    	    


class Spectrometer():
    "Spectrometer properties and methods container class " 
    def __init__(self, SerialNumber,tokamak):

	self.SerialNumber = SerialNumber
	self.Tokamak = tokamak
	

    #def fittingFunction( p,x): # ordinary gaussian
	#if len(p) != 3:
	    #raise Exception('WrongNumberOfParameters')
	
	#print 'spec fit fun'
	
	#p2 = empty(5)
	#p2[:3] = p
	#p2[3:] = 0
	#shape = UniversalLineShape( p2,x)
	#return shape
	
    #def lineCharacteristics(p):  # do nothing
	#return p
	
	
		
    def elementIdentification(self,wavelength,sigma):
	"""
	Find the closest line to the "wavelength" from lines in tokamak database. If any line was found, amed "unknown" is returned. 
	"""
	dispersion = (self.wavelength_range[1]-self.wavelength_range[0])/self.resolution
	presicion = min(sqrt(self.wavelength_presicion**2+sigma**2),dispersion*3)
	close_line_list = list()
	closest_line = inf
	element_name  = ""

	for i in self.Tokamak.elements:
	    index = abs(i[0]-wavelength) < presicion
	    if any(index):
		close_line_list.append([i[0][index], i[1]])
		print '\t\t\t',i[1],i[0][index], (i[0]-wavelength)[index]
	    close_line_i = argmin(  abs(i[0]-wavelength))
	    close_line = i[0][close_line_i]
	    
	    if abs(closest_line-wavelength)>abs(close_line-wavelength) and abs(close_line-wavelength) < presicion:
		closest_line = close_line
		element_name = i[1]	
		
	if isnan(closest_line):
	    print('\t\t\tclosest line: %4.3f , difference: %2.3f' % (wavelength, sigma))
	if closest_line == inf:
	    closest_line = wavelength
	    element_name = 'unknown'
	    
	    
	return element_name,closest_line,close_line_list
    
    def convertCountToPhotons(self,wavelength,intensity, make_plot = False):
	"""
	Estimate the number of photons necessery to achieve number of counts in "intensity" on the "wavelength". 
	This estimation is not based on real measurement of effectivity, but only as multiplying of effectivity of CCD,
	gratting and optical fiber. 
	"""
	f1 = interpolate.interp1d(self.sensor_efficiency[:,0],self.sensor_efficiency[:,1], bounds_error=True, fill_value=0.3) 
	f2 = interpolate.interp1d(self.gratting_efficiency[:,0],self.gratting_efficiency[:,1], bounds_error=True, fill_value=0.1) 
	fiber_effectivity = 10**(-self.fiber_attenuation[:,1]/10*self.fiber_length)
	f3 = interpolate.interp1d(self.fiber_attenuation[:,0],fiber_effectivity, bounds_error=True, fill_value=0.1) 
	#I dont't know sensitivity of CCD in UV !!
	if make_plot:	
	    plot(wavelength,self.photonsPerCount/(f1(wavelength)*f2(wavelength*f3(wavelength))))
	    xlabel('wavelength [nm]')
	    ylabel('photons/count')
	    show()
	    
	return intensity/(f1(wavelength)*f2(wavelength)*f3(wavelength))*self.photonsPerCount
    
    #def estimateNoiseInCounts(self,wavelength, intensity):
	#return ones(len(intensity))*self.readoutNoiseRMS
    #TODO ještě přetížit ty funkce na extrakci
    
    #TODO rovnou při extrakci určit readoutNoiseRMS
	
	

	
    def peakFitting(self,x0,ix0,wavelength,intensity, plot_fit = False):
	
	"""
	Identify the peak around the middle in index ix0 (wavelength x0). The peak is fitted by NLLSQ and centre of mass of line is compared with lines from database. Option plot_fit will show the plot with data, their best fit and spectral lines near of this line. 
	"""
	
	yn = intensity[ix0]
	ixr = ix0+1 
	n = len(wavelength)
	
	
	while ixr < n and yn > intensity[ixr]-sqrt(2):
	    yn = intensity[ixr]
	    ixr+=1
	yn = intensity[ix0]
	ixl = ix0
	while ixl > 0 and yn > intensity[ixl-1]-sqrt(2):
	    yn = intensity[ixl-1]
	    ixl-=1
	
	if ixr-ixl-3 <= 0:		#probably false peak, too much narrow
	    ixr = min(n, ixr+3)
	    ixl = max(0, ixl-3)


	a = intensity[ix0]*sqrt(2*pi*mean(self.inst_width[:,1])**2)
	p0 = (a,mean(self.inst_width[:,1]),x0)
	
	fun = lambda p,x: intensity[ixl:ixr]-self.fittingFunction(p,x)
	popt,pcov,info,mesg,success = leastsq(fun, p0, args = wavelength[ixl:ixr], full_output=1)
	
	if any(pcov == None):
	    pcov = ones((3,3))*nan

	# calculate real paramethers of this line
	popt_corr = self.lineCharacteristics(popt)
	
	
	chisq=sum(info["fvec"]*info["fvec"])
	doF=ixr-ixl-3
	print('\t\twavelength %4.3f +/- %1.3f; chi^2/doF: %4.2f'%(popt_corr[2],sqrt(pcov[2,2]),chisq/max(doF,1)))
	if chisq > scipy.stats.chi2.isf(0.1,doF):   # estimate uncorrect fit
	    pcov*=chisq/doF
	    
	# find closest lines
	line_name, line ,close_line_list= self.elementIdentification(popt_corr[2],sqrt(pcov[2,2]))

	#if   line == inf:
	    #plot_fit = True

	
	#plot graph for easier debugging
	if plot_fit:
	    step(wavelength[max(ixl-3,0):min(ixr+3,n)],intensity[max(ixl-3,0):min(ixr+3,n)],'r',where='mid')
	    step(wavelength[ixl:ixr+1]-(wavelength[ixl+1]-wavelength[ixl])/2, intensity[ixl:ixr+1],'b',where='post')  # steps over the integrating area of the pixel
	    interv = arange(max(ixl-3,0),min(ixr+3,n))
	    errorbar(wavelength[interv],intensity[interv], yerr=ones(len(interv)),fmt='r.')   #intensity with errorbars
	    errorbar(wavelength[ixl:ixr], intensity[ixl:ixr], yerr=ones(ixr-ixl),fmt='b.')  #intensity used for fitting with errorbars

	    interv   = linspace(wavelength[max(ixl-3,0)],wavelength[min(ixr+2,n-1)],100*(ixr-ixl))
	    plot(interv,self.fittingFunction(popt,interv), '--')

	    axhline(y=4, ls= '--',label = 'threshold')
	    
	    print('wavelength  %4.3f +/- %1.3f' % (popt_corr[2], sqrt(pcov[2,2])))
	    axvline(x = popt_corr[2], ls = '-.',c = 'g', label = "Center of mass") 

	    colour = ( 'g','r','c', 'm', 'y', 'k' )

	    for i, element in enumerate(close_line_list):
		for j,line in enumerate(element[0]):
		    axvline(x = line , label = element[1], c = colour[(i%len(colour))])
		    
	    axis('tight')
	    xlabel('$\lambda$ [nm]')
	    ylabel('SNR[-]')

	    legend(loc = 'best')
	    title('wavelength %4.3f +/- %1.3f; chi^2/doF: %4.2f'%(popt_corr[2],sqrt(pcov[2,2]),chisq/doF))
	    show()
	    close()
	    
	    
	d = diff(wavelength)[max(0,ix0-1)]

	ixl = max(int(ix0-popt_corr[1]/d*3), ixl)
	ixr = min(int(ix0+popt_corr[1]/d*3), ixr)
	return popt_corr,sqrt(diag(pcov)),ixl,ixr,line_name, line

    	
	
    def identifyLines(self,wavelength, intensity,readoutNoiseRMS,  debug = True):
	"""
	Detect line in the spectra given by field "intensity". Detection is based on chi2 statistics criterium. The detection threshold is set so that false detect of line in the spectrum with purely gaussian noise with sigma equal to readoutNoiseRMS is 1%. This method is more sensitive to weak lines than the naive detection of pixels higher than some threshold. 
	Lines are detected in infinite cycle until any line is detected.  The peak is fitted by function "peakFitting"  and at the and of the cycle the peaks are removed. 
	There are two ways how to remove peak. The peak can be substracted from the spectrum or the peak can be replaced by gaussian noise. Substrating can be used only of the instumental function is measured very precisely and therefore chi2 of the fit is close to 1. Because this is not fulfiled, tzhe second way is now used. 
	
	"""
	print '\t identifing lines'

	
	
	spectraParametres = list()
	line_fit = list()
	length = len(intensity)

	derivation = zeros(length)	
	
	mean_peak_FWHM = 6 #[px]  #TODO je to nesystematické, dát do vlastností spektrometru
    
	lim = scipy.stats.chi2.isf(0.01/length,mean_peak_FWHM/2)*2 #using only positive part (50% of data), misscorrect detect probability 1%

	#intensity-= median(intensity[self.black_pixels])   #  dávat to sem? 
	#readoutNoiseRMS = MAD(intensity)
	intensity/=readoutNoiseRMS	
	
	step = 0
	while True:   
	    step+= 1
	    if step >200:
		print('infinitely cycle')
		break
		
	    
	    cumSumSpectra = cumsum((intensity*(intensity>0))**2)*2 - arange(length)  #using only positive part (50% of data)
	    
	    derivation[2:length-2] = cumSumSpectra[4:]-cumSumSpectra[:length-4]+mean_peak_FWHM
	    
	    if all(derivation < lim):
		break


	    

	    #find maximum on the neighbourhood
	    ix0 = argmax(derivation)
	    ix0 = argmax(intensity[max(0,ix0-mean_peak_FWHM):min(ix0+mean_peak_FWHM, length)])+max(0,ix0-mean_peak_FWHM)
	    x0 = wavelength[ix0]
	    
	    ##if ix0 == self.resolution-1:
	    #plot(derivation, label= 'derivation')
	    ##plot((0,length), (lim,lim), label = 'threshold')
	    #axhline( y = lim, c = 'r', ls = '--',label = 'threshold')
	    #axvline( x = ix0, c = 'g',label = 'max')
	    #axvline( x = argmax(derivation), c = 'b', ls = '--',label = 'max deriv')
	    ##yscale('log', nonposy='clip')	
	    #axis('tight') 
	    #legend()
	    #show()
	
	    #plot(intensity, label= 'intensity')
	    #axvline( x = ix0, c = 'g',label = 'max')
	    #axvline( x = argmax(derivation), c = 'b', ls = '--',label = 'max deriv')
	    #axis('tight') 
	    #legend()
	    #show()
	
	    
		

	
	    if debug:
		plot(wavelength, cumSumSpectra, label = 'culumative sum of squars')
		axvline( x = x0, c = 'r', ls = '--')
		title('cumulative sum of squars of data')
		show()
		
	    
	    popt,sig,ixl,ixr,line_name,line = self.peakFitting(x0,ix0,wavelength,intensity,debug)
	    #renormalize
	    popt[0] = self.convertCountToPhotons(x0,readoutNoiseRMS*popt[0] )
	    sig[0] = self.convertCountToPhotons(x0,readoutNoiseRMS*sig[0] )
	    sig*=stats.t.isf(0.16,ixr- ixl-3) 	#"1 sigma" probability (68%)
	    
	    # wavelength, name, area,error,FWHM,error, Delta lambda, error
	    line_fit.append((line,line_name, popt[0],sig[0], popt[1],sig[1],(popt[2]-line ), sig[2]))

	    savetxt('current_peak.txt',(wavelength[ixl:ixr],intensity[ixl:ixr]))
	    # remove the fitted peak (in ideal case the peak should be substracted
	    #plot(wavelength[ixl+1:ixr-1], intensity[ixl+1:ixr-1],'o')
	    #plot(wavelength[ixl+1-5:ixr-1+5], intensity[ixl+1-5:ixr+5-1],'.r')
	    #show()
	    peakInterv = range(max(ixl+1,0),min(ixr, self.resolution-1)-1)
	    intensity[peakInterv] = random.randn(len(peakInterv))

	return line_fit
	
	
	
	
    
    
    
class OceanOptics_Spec(Spectrometer):
    "child of class Spectrometer with the Ocean Optics specific functions"
    def __init__(self, SerialNumber,tokamak):
	    Spectrometer.__init__(self, SerialNumber,tokamak)
	    #self.SerialNumber
	    #self.Tokamak
	    #self.calibrationPolynomeCorrection = zeros((4,1))
	    #self.wavelength_presicion =0.1#nm
	    #self.black_pixels = nan
	    #self.resolution = 3648
	    #self.min_integ_time = 3.8
	    #self.detector = 'TCD1304AP'
	    #self.gratting_efficiency = array(((0,1),(1200,1)))
	    
	    #self.sensor_efficiency = loadtxt('./HR4000/citlivost')
	    #self.photonsPerCount = 60
	    #self.hotPixels = None
	    
	    #self.inst_width = array(((0,0),(1200,0)))
	    #self.inst_skewness = array(((0,0),(1200,0)))
	    #self.voitig_parameter = 0
	    
	    
	
	    self.NonLinCoeff = zeros(7)		# implementovat to vůbec??
	    self.NonLinCoeff[0] = 1 
	    database = os.listdir('SpectralLines')


	    #load everything from config file
	    self.loadConfig()
	    
	    

    def start(self):
	"""
	It will start the java constrol program of the spectrometer. 
	"""


	def isThisRunning( process_name ):
	    ps     = subprocess.Popen("ps -eaf", shell=True, stdout=subprocess.PIPE)
	    output = ps.stdout.read()
	    output = output.find(process_name)
	    ps.stdout.close()
	    ps.wait()
	    
	    if output == -1:
		return False
	    else:
		return True
	
	java_library = '/opt/OmniDriver-/OOI_HOME/HighResTiming.jar:/opt/OmniDriver-/OOI_HOME/OmniDriver.jar'

	
	if isThisRunning('HighSpeedAcquisitionSample'):   #TODO dodělat
	    print 'spectrometer is already running'
	    return

	    
	try:
	    os.remove('nohup_java.out')
	    os.remove('highspeedacquisitionsample/HighSpeedAcquisitionSample.class')
	except:
	    print 'file can not be removed' 
	
	
	
	
	os.system('javac -classpath '+java_library+' HighSpeedAcquisitionSample.java')
	
	shutil.copyfile('HighSpeedAcquisitionSample.class', 'highspeedacquisitionsample/HighSpeedAcquisitionSample.class')
	# need admin rights to start with so high priority
	run_options = (self.SerialNumber,int(self.Tokamak.dataAcqusitionTime*1000), int(1000*self.integTime))
	os.system('nohup chrt -f 99 java -Djava.library.path=. -classpath '+java_library+':. highspeedacquisitionsample.HighSpeedAcquisitionSample %s %i %i & > nohup_java.out'  %run_options) 
	#os.system('tail -f nohup.out &')
	
	
	
	
		
	
	
    def wavelengthCorrection(self, old_wavelength):
	"""
	Apply the wavelength correction polynome
	"""
	new_wavelegth = copy(old_wavelength)
	for i,a in enumerate(self.calibrationPolynomeCorrection):
	    new_wavelegth+= a*old_wavelength**i
	
	return  new_wavelegth
	
	
	
	
    def loadConfig(self):
	"""
	Load the configuration file of the spectrometer
	"""
	def parseStringOfField(s):
	    s = s.splitlines() 
	    field = list()
	    for line in s:
		row = (line.strip('[]')).split()
		field.append(list())    
		for col in row:
		    field[-1].append(float(col))
	    return array(field)
	    
	    
	if not os.path.isfile(self.SerialNumber+'.cfg'):
	    raise Exception('Spectrometer config file does not exists.')	  
	    
	config = ConfigParser.RawConfigParser()	
	config.read(self.SerialNumber+'.cfg')
	
	
	
	SerialNumber = config.get('Basic parameters', 'Serial Number')
	
	if SerialNumber != self.SerialNumber:
	    raise Exception('WrongSerialNumber')
	
	self.detector = config.get('Basic parameters', 'Detector') 
	self.GratingTyp = config.get('Basic parameters', 'Optical grating') 
	self.wavelength_presicion = config.getfloat('Basic parameters', 'Wavelength presicion [nm]') 
	self.min_integ_time = config.getfloat('Basic parameters', 'Minimum integration time [ms]') 
    	tmpstr = config.get('Basic parameters', 'Opticaly black pixels')    
	self.black_pixels = int_(parseStringOfField(tmpstr))
	self.resolution = config.getint('Basic parameters', 'Resolution [px]' )    
	self.photonsPerCount = config.getfloat('Basic parameters', 'min photons/count')  
	
	tmpstr = config.get('Efficiencies', 'Gratting efficiency')
	self.gratting_efficiency = parseStringOfField(tmpstr)
	tmpstr = config.get('Efficiencies', 'Sensor efficiency')  # pro HR4000 to velmi výrazně nesedí s údaji na webu
	self.sensor_efficiency = parseStringOfField(tmpstr)
	tmpstr = config.get('Efficiencies', 'optical fibers attenuation (dB/m)')  
	self.fiber_attenuation = parseStringOfField(tmpstr)	
	self.fiber_length = config.getfloat('Efficiencies', 'fiber length') 
	
	tmpstr = config.get('Other properties', 'Calibration polynome correction')
	self.calibrationPolynomeCorrection = parseStringOfField(tmpstr)[0]
	tmpstr = config.get('Other properties', 'Linearity polynome correction')
	self.linearityPolynomeCorrection = parseStringOfField(tmpstr)[0]
	
	tmpstr = config.get('Other properties', 'Hot pixels')
	self.hotPixels = parseStringOfField(tmpstr)[0]	
	tmpstr = config.get('Other properties', 'wavelength range') 
	self.wavelength_range = parseStringOfField(tmpstr)[0]	
	
	tmpstr = config.get('Instrumental broadening function', 'width(wavelength)') 
	self.inst_width = parseStringOfField(tmpstr)
	self.F_inst_width = interpolate.interp1d(self.inst_width[:,0],self.inst_width[:,1],bounds_error=False,fill_value=0) 
	tmpstr = config.get('Instrumental broadening function', 'skewness(wavelength)') 
	self.inst_skewness = parseStringOfField(tmpstr)
	self.F_inst_skewness = interpolate.interp1d(self.inst_skewness[:,0],self.inst_skewness[:,1],bounds_error=False,fill_value=0) 
	self.voitig_parameter = config.getfloat('Instrumental broadening function', 'voitig parameter') 

	
	
	print 'Loaded data of ', SerialNumber

	
    def saveConfig(self):
	"""
	Save the configuration file of the spectrometer
	"""    
	config = ConfigParser.RawConfigParser()
	
	config.add_section('Basic parameters')
	config.set('Basic parameters', 'Type', 'Ocean Optics Spectrometer')
	config.set('Basic parameters', 'Serial Number', self.SerialNumber )
	config.set('Basic parameters', 'Detector', self.detector )
	config.set('Basic parameters', 'Optical grating', self.GratingTyp )

	config.set('Basic parameters', 'Date', asctime())	
	config.set('Basic parameters', 'Wavelength presicion [nm]', self.wavelength_presicion) 
	config.set('Basic parameters', 'Opticaly black pixels', self.black_pixels )    
    	config.set('Basic parameters', 'Resolution [px]', self.resolution )    
	config.set('Basic parameters', 'Minimum integration time [ms]', self.min_integ_time)  
	config.set('Basic parameters', 'min photons/count', self.photonsPerCount)  
	config.add_section('Efficiencies')
	config.set('Efficiencies', 'Gratting efficiency', self.gratting_efficiency )	
	config.set('Efficiencies', 'Sensor efficiency', self.sensor_efficiency)
	config.set('Efficiencies', 'Optical fibers attenuation (dB/m)', self.fiber_attenuation)
	config.set('Efficiencies', 'fiber length', self.fiber_length)
	
	
	config.add_section('Other properties')
	config.set('Other properties', 'Calibration polynome correction', self.calibrationPolynomeCorrection)
	config.set('Other properties', 'linearity polynome correction'  , self.linearityPolynomeCorrection)
	config.set('Other properties', 'Hot pixels', self.hotPixels)
	config.set('Other properties', 'wavelength range',self.wavelength_range ) 

	config.add_section('Instrumental broadening function')
	config.set('Instrumental broadening function', 'width(wavelength)', self.inst_width )  # šířka kalibračních čar v několika bodech	
	config.set('Instrumental broadening function', 'skewness(wavelength)', self.inst_skewness )  # šikmost kalibračních čar v několika bodech	
	config.set('Instrumental broadening function', 'voitig parameter', self.voitig_parameter )  # negausovkost základny	

	with open(self.SerialNumber+'.cfg', 'wb') as configfile:
	    config.write(configfile)
    
    def checkCalibration(self):  # vykresí to graf a určí chi^2  = suma(delta labda ku sigma)^2
	#zapolá to makeCalibration
	pass
    
    
    def makeCalibration(self, dryRun = True):
	pass
	#načte list kalibračních čar
	#100x spustí sběr data a zprůměruje
	#projde kalibrační čáry a najde nejbližší peak
	#aktualizuje (pro dryRun = false) to self.inst_width, self.inst_skewness, calibrationPolynomeCorrection, 
	# vykreslí to současné a nové hodnoty (i errorbary)
	#TODO nastavit width(wavelength),skewness, voitig parameter, Calibration polynome correction
	
    def calcParamsForFitting(self,p):
	"""
	The lines are not fitted by full set of paramathes of the general peak function UniversalShape. Only a few paramethers are used, other are calculated here. 
	""" 
	p_full = empty(5)
	p_full[:3] = p[:3]
	#print 'int_w ', self.F_inst_width(p[2]), self.F_inst_skewness(p[2]), self.voitig_parameter 
	#p_full[1] = p[1]
	p_full[1] = sqrt(p[1]**2+self.F_inst_width(p[2])**2)  #add instrumental broadening

	p_full[3] = self.F_inst_skewness(p[2])
	#p_full[3] = 10
	p_full[4] = self.voitig_parameter 
	#p_full[4] = 0.00
	return p_full
	
    def fittingFunction(self,p,x):  # využívat informace o znalosti instrumentálního rozšíření
	p_full = self.calcParamsForFitting(p)	
	return UniversalLineShape( p_full,x)

    def lineCharacteristics(self,p):
	"""
	The lines are not fitted by full set of paramathes of the general peak function UniversalShape. 
	Therefore the real properties of the peak are calcutelad here from the limited set of the paramethers. 
	
	""" 
	
	p_full = self.calcParamsForFitting(p)
	p_corrected = lineCharacteristics(p_full)
	p_corrected[1] = sqrt(p_corrected[1] **2-self.F_inst_width(p_corrected[2])**2)
	return lineCharacteristics(p_full)
    

	
	
    def loadData(self,path,typ,name = 'spectra'):
	"""
	Load data from file and create object "shot" containig all necessary informations andou the tokamak shot. It supports many different types of files. 
	Java_Data_file  - file experted from java control program HighSpeedAcquisitionSample.java
	Python_Data_file - file exported from this python program
	SpectraSuite_xml - file exported from OceanOptics SpectralSuite
	SpectraSuite_txt - acsii file exported from OceanOptics SpectralSuite
	SpectraSuite_HighSpeed - ascii file exported from OceanOptics SpectralSuite after high speed acquisition
	"""
	
	shot = Shot(self.Tokamak.shotNumber,  self)
	
	if typ == 'Java_Data_file':   # acqurired by HighSpeedAcquisitionSample.java
	    #extract data about spectra

	    
	    shot.DataFile = path+name+'.txt'

	    if not os.path.isfile(shot.DataFile):
		raise Exception('Data file '+shot.DataFile+' does not exists.')
	    
	    
	    f = open(shot.DataFile, 'r')

	    for i,line in enumerate(f):
		if line[:14] == 'Date and time:':
		    shot.time = mktime(strptime(line[15:-1], "%Y.%m.%d %H:%M:%S"))
		if line[:18] == 'Number of spectra:':
		    shot.n_spectra = int(line[20:-1])
		if line[:25] == 'Number pixels in spectra:':
		    shot.n_pixels = int(line[27:-1])
		if line[:22] == 'Integration time [us]:':
		    shot.integ_time = int(line[24:-1])/1000.0
		if line[:22] == 'Board temperature [C]:':
		    shot.temperature = float(line[24:-1])
		if line[:17] == 'Spectrometer S/N:':  
		    SN = line[18:-1]
		if line[:18] == 'Exact time stamps:':
		    break
	    
	    if SN != self.SerialNumber:
		raise Exception('wrong serial number, '+SN+' != '+self.SerialNumber)

		
	    shot.time_stamps = zeros(shot.n_spectra)
	    for i,line in enumerate(f):
		if i >= shot.n_spectra:
		    break
		m = re.search('\s[0-9]+\.[0-9E]+',line)
		s = m.group(0)
		shot.time_stamps[i] = float(s)
	    shot.time_stamps/= 1e6   #convert ns -> ms	
	    shot.time_stamps+= -shot.time_stamps[0]+self.Tokamak.trigger_shift #shift ->0	    
	    
	    data = genfromtxt(f)
	    f.close()

	    shot.wavelength = data[:shot.n_pixels]
	    shot.readoutNoiseRMS = 0
	    for i in range(shot.n_spectra):
		intensity = data[(i+1)*shot.n_pixels: (i+2)*shot.n_pixels]
		#intensity -=  median(intensity[self.black_pixels])
		shot.spectra.append(Spectrum(shot,shot.wavelength,intensity,shot.time_stamps[i] ))
		shot.readoutNoiseRMS  += std(intensity[self.black_pixels])**2
	    shot.readoutNoiseRMS = sqrt(shot.readoutNoiseRMS/shot.n_spectra)  # Zkontrolovat stabilitu odhadu
	    #extract data about background and dark current	
	    if  os.path.isfile(path+'background.txt'):	
		imax = 0
		f = open(path+'background.txt', 'r')
		for i,line in enumerate(f):
		    if line[:-1].replace('.','').isalnum():
			imax = i
			break
		    if line[0:22] == 'Integration time [us]:':
			integ_time_background = int(line[24:-1])/1000
		f.close()
	
		shot.darkNoise = genfromtxt(path+'background.txt', skip_header=imax)
		shot.darkNoise -= median(shot.darkNoise[self.black_pixels])	
		shot.darkNoise*=shot.integ_time/float(integ_time_background)
	    

	    
	    #extract data about readout patterns
	    if  os.path.isfile(path+'readoutpatterns.txt'):
		imax = 0
		f = open(path+'readoutpatterns.txt', 'r')
		for i,line in enumerate(f):
		    if line[:-1].replace('.','').isalnum():
			imax = i
			break
		f.close()
		
		shot.readOutPatterns =  genfromtxt(path+'readoutpatterns.txt', skip_header=imax)
		shot.readOutPatterns-= (shot.darkNoise/shot.integ_time)*self.min_integ_time
	    else:
		shot.readOutPatterns = mean(intensity[self.black_pixels])    
		
		
	if typ == 'Python_Data_file':
	    	    
	    shot.DataFile = path+name+'.txt'

	    if not os.path.isfile(shot.DataFile):
		raise Exception('Data file '+shot.DataFile+' does not exists.')
	    
	      
	    f = open(shot.DataFile, 'r')
	    
	    for i,line in enumerate(f):
		ind = line.find('Serial Number')
		if ind != -1:
		    SN =  line[ind+15:-1]
		ind = line.find('Date and time (GMT)')
		if ind != -1:    
		    shot.time =  mktime(strptime(line[ind+21:-1]))
		ind = line.find('Number of spectra')
		if ind != -1:
		    shot.n_spectra = int(line[ind+19:-1])
		ind = line.find('Resolution')
		if ind != -1:
		    shot.n_pixels = int( line[ind+16:-1])  
		ind = line.find('Integration time [ms]')
		if ind != -1:
		    shot.integ_time = float(line[ind+22:-1])
		ind = line.find('Board temperature [C]')
		if ind != -1:
		    shot.temperature = float(line[ind+22:-1])
		ind = line.find('Time stamps [ms]')
		if ind != -1:
		    s = line[ind+18:-1]
		    shot.time_stamps = float_(s.strip('[]').split())
		ind = line.find('Noise RMS')
		if ind != -1:	
		    print line[ind+20:-1]
		    shot.readoutNoiseRMS = float(line[ind+20:-1])
		if line.find('***************') != -1:
		    break
	    if SN != self.SerialNumber:
		raise Exception('wrong serial number, '+SN+' != '+self.SerialNumber)

		
	    data = genfromtxt(f)
	    f.close()
	    shot.wavelength = data[:,0]
	    
	    for i in range(shot.n_spectra):
		shot.spectra.append(Spectrum(shot,shot.wavelength, data[:,i+1], shot.time_stamps[i]))
	    
	    shot.darkNoise = zeros(shot.n_pixels)
	    shot.readOutPatterns = zeros(shot.n_pixels)

	    
	    
	if typ == 'SpectraSuite_xml':		#acqurid by Data Studio

	    shot.temperature = nan
	    shot.DataFile = path+name+'.ProcSpec'

	    if not os.path.isfile(shot.DataFile):
		raise Exception('Data file '+shot.DataFile+' does not exists.')
	    
	    
	    if  os.path.isfile('./readoutpatterns.txt'):
		shot.readOutPatterns =  genfromtxt('./readoutpatterns.txt', skip_header=1)
		shot.readOutPatterns -= median(shot.readOutPatterns[self.black_pixels])   
	    
	    
	    with zipfile.ZipFile(shot.DataFile, 'r') as myzip:
		shot.n_spectra = len(myzip.namelist())-2

		for data in myzip.namelist():

		    if data !=  'OOISignatures.xml' and data != 'OOIVersion.txt':
			#print directory+'/'+data_file+'/'+data
			
			t = myzip.getinfo(data).date_time
			shot.time = mktime((t[0],t[1],t[2],t[3],t[4],t[5],0, 0,0))

			
			lines = myzip.open(data)
			
			
			for j,line in enumerate(lines):
			    if line.find('<numberOfPixels>') != -1:
				ind = line.find('<numberOfPixels>')
				shot.n_pixels =  int(line[ind+16:-18])			

			    if line.find('<integrationTime>') != -1:
				ind = line.find('<integrationTime>')				
				shot.integ_time =  int(line[ind+17:-19])/1000.0#convert to ms
			    if line.find('<spectrometerSerialNumber>') != -1:
				ind = line.find('<spectrometerSerialNumber>')
				SN  = line[26+ind:-28]				
				

			if SN != self.SerialNumber:
			    raise Exception('wrong serial number, '+SN+' != '+self.SerialNumber)
			
			lines.close()
			
			lines = myzip.open(data)
			j = 0

			intensity = zeros(shot.n_pixels)
			shot.wavelength  = zeros(shot.n_pixels)
			channelWavelengths = False
			sourceSpectra = False
			pixelValues = False

			for i,line in enumerate(lines):
			    index = line.find('<double>')
			    
			    if line.find('<com.oceanoptics.omnidriver.spectra.OmniSpectrum>') != -1:
				sourceSpectra = True
				j = 0	
				
			    if line.find('</com.oceanoptics.omnidriver.spectra.OmniSpectrum>') != -1:
				sourceSpectra = False	
				
			    if line.find('<channelWavelengths>') != -1:
				channelWavelengths = True
				j = 0    
				
			    if line.find('</channelWavelengths>') != -1:
				channelWavelengths = False    
				
			    if line.find('<pixelValues>') != -1:
				pixelValues = True
				j = 0    
				
			    if line.find('</pixelValues>') != -1:
				pixelValues = False      
			    
				
			    if index != -1:
				if sourceSpectra and pixelValues :
				    intensity[j] =  float(line[index+8:-10])
				if sourceSpectra and channelWavelengths:
				    shot.wavelength[j] = float(line[index+8:-10])
				j+=1
			lines.close()
			
			intensity -=  median(intensity)

			shot.spectra.append(Spectrum(shot,shot.wavelength,intensity)) 

	       
	    
	    shot.time_stamps = arange(shot.n_spectra)*shot.integ_time+self.Tokamak.trigger_shift 
	    shot.readoutNoiseRMS = MAD(intensity)
	   
	   
	if typ == 'SpectraSuite_txt':		#text acqurid by Data Studio
	    shot.temperature = nan
	    shot.DataFile = path+name+'.txt'
	    if not os.path.isfile(shot.DataFile):
		raise Exception('Data file '+shot.DataFile+' does not exists.')
	    
	    
	    if  os.path.isfile('./readoutpatterns.txt'):
		shot.readOutPatterns =  genfromtxt('./readoutpatterns.txt', skip_header=1)
		shot.readOutPatterns -= median(shot.readOutPatterns[self.black_pixels])   
		
	    f = open(shot.DataFile, 'r')
	    line = f.readline()
	    if line != 'SpectraSuite Data File\n':
		raise Exception('wrong file type')
	    
	    while True:
		line = f.readline()
		if not line:
		    if not line: raise Exception('corrupted header of file')
		if line[:5] == 'Date:':
		    try:
			shot.time = mktime(strptime(line[6:-1], "%a %b %d %H:%M:%S %Z %Y")) # sometimes it can make mysterious errors connected with figure()
		    except:
			shot.time = os.path.getctime(shot.DataFile)
		if line[:36] == 'Number of Sampled Component Spectra:':
		    shot.n_spectra = int(line[37:-1])
		if line[:16] == 'Number of Pixels':
		    shot.n_pixels = int(line[-5:-1])
		if line[:24] == 'Integration Time (usec):':
		    if line[-11:-1] == '('+self.SerialNumber+')':
			shot.integ_time = int(line[25:-11])/1000.0 #convert to ms	
		    else:
			shot.integ_time = int(line[25:-1])/1000.0 #convert to ms	
		if line[:14] == 'Spectrometers:' or line[:27] == 'Spectrometer Serial Number:':  #TODO opravit od do
		    SN = line[-9:-1]
		    if SN != self.SerialNumber:
			raise Exception('wrong serial number, '+SN+' != '+self.SerialNumber)
		    
		if line[:5] == '>>>>>':
		    break
	    
	    if shot.n_spectra == 0:  # missing information in file
		shot.n_spectra = 1
	    


	    
	    replaceComma = lambda s:float(s.replace(',','.'))
	    data = genfromtxt(f, skip_footer = 1, converters = {0:replaceComma, 1:replaceComma })
	    f.close()
	    shot.wavelength = data[:,0]
	    intensity = data[:,1]
	    intensity -=  median(intensity)
	    shot.spectra.append(Spectrum(shot,shot.wavelength,intensity)) 
	    shot.readoutNoiseRMS  = MAD(intensity)
	    #shot.readoutNoiseRMS  = std(intensity[self.black_pixels])

	   
	if typ == 'SpectraSuite_HighSpeed':   
	
	    
	    shot.temperature = nan
	    shot.DataFile = path+name+'.txt'
	    if not os.path.isfile(shot.DataFile):
		raise Exception('Data file '+shot.DataFile+' does not exists.')
	    
	    f = open(shot.DataFile, 'r')
	    
	    first_line = f.readline()
	    if first_line == 'SpectraSuite Data File\n':
		raise Exception('wrong file type')
	    
	    data = loadtxt(f)
	    f.close()
	    
	    shot.time_stamps = float_(first_line.split())+self.Tokamak.trigger_shift 
	    shot.time = os.path.getctime(shot.DataFile) # it can be wrong time!
	    shot.integ_time = shot.time_stamps[-1]/(len(shot.time_stamps)-1)  # only estimation!
	    shot.n_spectra = size(data,1)-1
	    shot.n_pixels = size(data,0)	    
	    shot.wavelength = data[:,0]
	    
	    for i in range(shot.n_spectra):
		intensity = data[:,i+1]
		intensity -= median(intensity)
		shot.spectra.append(Spectrum(shot,shot.wavelength,intensity,shot.time_stamps[i] )) 
	    

	    shot.readoutNoiseRMS  = MAD(data[self.black_pixels,1:])

	    if max(abs(amin(shot.wavelength)-self.wavelength_range[0]), abs(amax(shot.wavelength)-self.wavelength_range[1]))>10:
		raise Exception('bad wavelength range %d - %d != %d - %d' %(amin(shot.wavelength),amax(shot.wavelength),self.wavelength_range[0],self.wavelength_range[1]))
	
	
	    
	if shot.n_pixels != self.resolution:
	    raise Exception('wrong number of pixels, '+str(shot.n_pixels)+' != '+str(self.resolution))
	    
	    
	return shot
	
	
	    
	
	
##TODO ve fitovací funkci zahrout i to , že je to přenormované šumem. (prostě vypočítávat šum i pro tu teoretickou funkci)    
#class CIII_HR_Spec(Spectrometer):   
    ##TODO různé konfigy prop různá nastavení? 
    #def __init__(self):
	#Spectrometer.__init__(self, SerialNumber,tokamak)
	##TODO musí si to pamatovat všechny nastavitelné parametry
	
    #def convertCountToPhotons(self,wavelength, intensity, make_plot = False):
	##TODO tady by si to mohlo pamatovat ten šum 
	##TODO založit to na vzorci z výzkumáku
	
	#if len(wavelength) != len(intensity):
	    #wavelength = ones(len(intensity))*mean(wavelength)
	    
	#if make_plot:
	    #plot(wavelength,f1(wavelength)*f2(wavelength)*self.photonsPerCount)
	    #xlabel('wavelength [nm]')
	    #ylabel('photons/count')
	    #show()    
	##TODO střeva
	#return corr_intensity # šum i intenzita přepočtené na fotony
	
    #def estimateNoiseInCounts(self,wavelength, intensity):
	#if len(wavelength) != len(intensity):
	    #wavelength = ones(len(intensity))*mean(wavelength)
	    
	    ##TODO dodělat střeva
	    
	#return noise  #šum v počtu countů
    #def PeakFitting:
    #def fittingFunction:# buď 
	## přidá navíc paramer základna, využije se tu maxumum dosavadních znalostí. 
    #def identifyLines(self,wavelength, intensity):# return table of lines, wavelngths, +other parametres


    
    
    
    
    
class  Shot():
    """Discharge container class

    used for plotting and obtaining discharge data
    """
    def __init__(self, shot_num,  Spectrometer):
	
	self.Spectrometer = Spectrometer
	self.shot_num=shot_num

	self.spectra = list()

	self.wavelength = nan  # udělat korekci při načítání
	self.readOutPatterns = zeros(Spectrometer.resolution) #odečíst mediaán!!!
	self.darkNoise = zeros(Spectrometer.resolution) #odečíst mediaán!!!, normované na integrační dobu
	self.readoutNoiseRMS = 13 # TODO co to je? chovám se k tomu jak o k poli
	self.maxSNR  = 0 #TODO co stím?
	self.n_spectra = 0
	self.time_stamps = zeros(1)
	self.DataFile = ''
	
   
    
    def plasmaShot(self): 
	"""
	Detect if there was any plasma during the shot. 
	"""
    
	plasma = zeros(self.n_spectra, dtype = 'bool')
	for i in arange(self.n_spectra):
	    plasma[i] = self.spectra[i].plasma()	
	return any(plasma)
	

    def getData(self,sensitivity_correction = False,noiseCorrected = True): 
	"""
	Return the field containig spectra from the all acqusitions during the shot.  Option sensitivity_correction set that the correction on the not constant sensitivity of the whole device at different wavelength shold by applied. Option noiseCorrected set that the readout noise patterns, dark noise and hot pixels should be removed. 
	"""
	spectra = empty((self.Spectrometer.resolution, self.n_spectra))
	
	for i in arange(self.n_spectra):
	    spectra[:,i] = self.spectra[i].getData(sensitivity_correction , noiseCorrected )
	
	return copy(self.wavelength), spectra
	
	
	
    def saveData(self,filename = None, noiseCorrected = True): 
	"""
	Export data from the shot in well arranged shape. At the bigining of the file is the header and than coloumns with wavelengths and intensities. 
	"""
	
	print 'saving data '#, (self.time_stamps[-1]/(len(self.time_stamps)-1))
	if filename == None:
	    filename = self.DataFile+'_Data_.txt'
	
	# make a header
	header = ''
	header += 'Spectrometer:\t'	+	'Ocean Optics Spectrometer\n'
	header += 'Serial Number:\t'	+ 	self.Spectrometer.SerialNumber+'\n'
	header += 'Date and time (GMT):\t'	+ 	asctime(localtime(self.time))+'\n'
	header += 'Number of spectra:\t'+	str(self.n_spectra)+'\n'
	header += 'Resolution [px]:\t'	+	str(self.Spectrometer.resolution)+'\n'
	header += 'Integration time [ms]:\t'+	'%2.2f' %(self.integ_time)+'\n'
	header += 'Acquisition time [ms]:\t'+	'%2.2f' %((self.time_stamps[-1]-self.time_stamps[0])/(len(self.time_stamps)-1))+'\n'
	header += 'Board temperature [C]:\t'+ 	'%2.3f' %self.temperature +'\n'
	header += 'Time stamps [ms]:\t'	+ 	str(self.time_stamps)+'\n'
	header += 'Noise RMS [counts]:\t'	+	'%2.1f' %self.readoutNoiseRMS+'\n'
	header += '*'*100+'\n'
	
	#save spectra
	spectra = empty((self.Spectrometer.resolution, self.n_spectra+1))
	spectra[:,0] = self.wavelength
	for i in arange(self.n_spectra):
	    spectra[:,i+1] = self.spectra[i].getData(noiseCorrected = noiseCorrected)
	    
	    
	    
	f = open(filename, 'w')
	f.write(header)
	savetxt(f, spectra,fmt=['%1.3f']+['%4i']* self.n_spectra)
	f.close()
    

    
	
	
	
    def plot(self,  file_type = 'pdf', log_scale = True ,sensitivity_correction = True, w_interval = None):
	"""
	Plot spectra in separate files. 
	Options:
	file_type  -  file type of exported graph
	log_scale - set log scale
	sensitivity_correction -  correction on the not constant sensitivity of the whole device at different wavelength 
	w_interval - interval if the plotted wavelengths
	"""
	
	print 'plotting separate graphs'
	
	for i,spectrum in enumerate(self.spectra):   	    
	    self.spectra[i].plot(sensitivity_correction, log_scale =log_scale, noiseCorrected = True, w_interval = w_interval)
	    if sensitivity_correction:
		ylabel('photons [-]')
	    else:
		ylabel('counts [-]')
		
	    if file_type == None:
		show()
	    else:
		savefig(self.DataFile+'_Graph_%d.'%(i)+file_type)
	    close()
	
	
	
    
    def plot_all(self, file_type = 'pdf', log_scale = True ,sensitivity_correction = True, w_interval = None, plasma = True):
	"""
	Plot spectra in to one plot
	Options:
	file_type  -  file type of exported graph
	log_scale - set log scale
	sensitivity_correction -  correction on the not constant sensitivity of the whole device at different wavelength 
	w_interval - interval if the plotted wavelengths
	plasma  - only spectra with some lines
	"""
	
	
	print 'plotting all together'
	    
	class MyFormatter(ScalarFormatter): 
	    def __call__(self, x, pos=None): 
		if pos==0: 
		    return '' 
		else: return ScalarFormatter.__call__(self, x, pos)    
	    

	
	fig = figure(num=None, figsize=(12, 12), dpi=80, facecolor='w', edgecolor='k')

	subplots_adjust(wspace=0,hspace=0)
	
	if not self.plasmaShot():
	    plasma = False
	    
	n_plots =  self.n_spectra
	
	for i,spectrum in enumerate(self.spectra):
	    if not spectrum.plasma() and plasma:
		n_plots-= 1
	
	i_plot = 0
	for i,spectrum in enumerate(self.spectra):
	    if (not spectrum.plasma()) and plasma:
		continue 
	    i_plot+= 1
	    ax = fig.add_subplot(n_plots,1,i_plot)
	    ax.xaxis.set_major_formatter( NullFormatter() )
	    ax.yaxis.set_major_formatter( MyFormatter() )
	    #subplot( self.n_spectra ,1, i+1)
	    
	    self.spectra[i].plot(sensitivity_correction, log_scale =log_scale, noiseCorrected = True, w_interval = w_interval)
	    if sensitivity_correction:
		ylabel('photons [-]', rotation='horizontal')
	    else:
		ylabel('counts [-]', rotation='horizontal')
	    
	    
	ax = fig.add_subplot(n_plots ,1,n_plots)  
	ax.xaxis.set_major_formatter( ScalarFormatter() )
	
	
	if file_type == None:
	    show()
	else:
	    savefig(self.DataFile+'_Graph_.'+file_type)
	close()

    def processShot(self, details = False):
	"""
	prepare table of identified lines and their properties for another functions
	"""
	
	
	LinesTable = []

	
	def addNewLine(wavelength,name):
	    LinesTable.append(list((wavelength,name)))
	    
	    if details: 
		emptyCell = list((0,nan,nan,nan,nan,nan))
	    else:
		emptyCell = list((0,))
		
	    for i in range(self.n_spectra):
		LinesTable[-1].append(emptyCell)
		
	def setValue(wavelength,time_slice, data):
	    for row in LinesTable:
		if row[0] ==  wavelength:
		    if details:
			row[2+time_slice] = data
		    else:
			row[2+time_slice] = data[0]
		    break
	    
	
	for i,spectrum in enumerate(self.spectra):   
	    Table = spectrum.processSpectrum()
	    print 'processing '+str(i)+'-th spectrum'

	    for line in Table:
		exist = False
		for row in LinesTable:
		    if row[0] ==  line[0]:
			exist = True
			setValue(line[0],i,line[2:] )
			break
		if not exist:
		    addNewLine(line[0],line[1])
		    setValue(line[0],i,line[2:] )

	ColoumnNames = list(('wavelength','line'))
	if details: 
	    OtherColoumns = list(('area','err','FWHM','err', 'Delta lambda','err'))
	else:
	    OtherColoumns = list(('area',))
	
	for i in range(self.n_spectra):
	    for name in OtherColoumns:
		ColoumnNames.append(name)
		
	return ColoumnNames,LinesTable,self.time_stamps
	    
	    
    def plotLinesTimeEvolution(self,file_type = None, only_identified = True, n_lines = 6):
	"""
	Identify annd plot the most intensive lines and plot their time evolution during the shot
	Options:
	file_type - file type of exported graph
	only_identified - extract only identified lines
	n_lines - number of plotted lines
	
	"""
	
	print 'lines evolution plotting'
	IntensitiesTable = []
	ErrorsTable = []
	NamesTable = []
	ColoumnNames,LinesTable,time_stamps = self.processShot(True) 
	
	
	for line in LinesTable:	    
	    if only_identified == True and line[1] == 'unknown':
		continue	    
	    NamesTable.append(line[1]+'  %3.1f' %line[0]+'nm')	    
    
	    intensity = zeros(self.n_spectra)
	    errorbars= zeros(self.n_spectra)
	    for j, time_slice in  enumerate(line[2:]):
		if not self.spectra[j].plasma():
		    continue
		intensity[j] = time_slice[0]
		errorbars[j] = time_slice[1]
	    IntensitiesTable.append(intensity)
	    ErrorsTable.append(errorbars)
	
	if len(IntensitiesTable) == 0:
	    return
	
	plasma_end = self.n_spectra
	for i in arange(self.n_spectra,0,-1):
	    if self.spectra[i-1].plasma():
		plasma_end = i
		break 
	
	IntensMax =  amax(IntensitiesTable, axis = 1)

	sort_ind = argsort(IntensMax, order=None)[::-1]	
	
	for i,index in enumerate(sort_ind):
	    if i > n_lines:
		break
	    errorbar(time_stamps,IntensitiesTable[index],yerr=ErrorsTable[index],label=NamesTable[index])

	legend(loc = 'best')
	t_start= time_stamps[0]
	t_end = time_stamps[min(plasma_end,self.n_spectra-1)]
	axis([t_start,t_end,0, amax(IntensitiesTable) *1.1])
	xlabel('t [ms]')
	ylabel('Intensity [a.u.]')
	title('The most intense lines')
	if file_type == None:
	    show()
	else:
	    savefig(self.DataFile+'_LinesEvolution_.'+file_type)
	close()

    def positionHistorogram(self):
	"""
	plot the hytorogram of positions of identified lines. It is useful for false line identifications. 
	"""
	
	ColoumnNames,LinesTable,time_stamps = self.processShot(True) 
	hist_list = list()
	for line in LinesTable:
	    for cell in line[2:]:
		if not isnan(cell[4]):
		    hist_list.append(( line[0] ,cell[4]))
	
	hist_list = array(hist_list)

	plot(hist_list[:,0], hist_list[:,1], '.')
	show()
	
	    
    def tableOfResultes(self,sufix = '',  sort = False, details = False, only_identified = True):
	"""
	Export table with lines and their properties
	"""
	
	ColoumnNames,LinesTable,time_stamps = self.processShot(details) 
	
	f = open(self.DataFile+'_Resultes_'+sufix+'.csv', 'w')
	
	
	
	for name in ColoumnNames:
	    f.write(name)
	    f.write('   ;')

	    
	    
	wavelengths = []	
	for i,row in enumerate(LinesTable):
	    wavelengths.append(row[0])
	  
	  
	if sort:
	    wavelengths.sort()
	
	  
	for pos in wavelengths:   # complicate way to sort the lines
	    
	    for row in LinesTable:
		if only_identified  and row[1] == 'unknown':
		    continue

		if row[0] == pos:
		    f.write('\n')
		    for cells in row:
			col = 0   
			if isinstance(cells,list) or isinstance(cells,tuple):
			    for num in cells:
				f.write(locale.format("%5.3f",num)+'   ;')
				col += 1

			else:
			    if isinstance(cells, double):
				f.write(locale.format('%3.3f', cells)+'   ;')
			    else:
				f.write(str(cells)+'   ;')
			    col += 1
		
		
	f.closed
	    
	    
	    
	    
	    
	    
	    
	    
	   	    

class  Spectrum():
    """
    Spectrum container class  used for plotting and preparing of single spectrum
    """
    
    def __init__(self,Shot,wavelength,intensity, time_stamp = 0):
	self.Shot = Shot
	self.intensity = intensity  #odečíst mediaán!!!, kontrolovat jeszli to neovlivňují další funkce
	self.wavelength  = Shot.Spectrometer.wavelengthCorrection(wavelength)
	
	self.time_stamp = time_stamp
	self.lineFitList = []
	
	
    def plasma(self):
	"""
	Detect if in the data are some lines    
	"""
	if any(self.getData() > self.Shot.readoutNoiseRMS*5):
	    return True
	return False
	
	
    def plot(self,sensitivity_correction = True, log_scale = True, noiseCorrected = True, w_interval = None):
	"""
	Plot a single spectrum
	log_scale - set log scale
	sensitivity_correction -  correction on the not constant sensitivity of the whole device at different wavelength 
	w_interval - interval if the plotted wavelengths
	"""
	if not w_interval:
	    w_interval = [amin(self.wavelength), amax(self.wavelength)]
	    
	    
	w_interval_bool = (self.wavelength >= w_interval[0]) *  (self.wavelength <= w_interval[1] )

	intensity = self.getData(sensitivity_correction,noiseCorrected)*w_interval_bool
	
	noise = self.Shot.readoutNoiseRMS	
	threshold = scipy.stats.norm.isf(0.05/self.Shot.Spectrometer.resolution)   # 5% probability of mistake
	if sensitivity_correction:
	    noise = self.Shot.Spectrometer.convertCountToPhotons(self.wavelength, noise)
	else:
	    noise *= ones(self.Shot.Spectrometer.resolution)
	
	min_threshold = amin(threshold*noise)
	noise *= w_interval_bool   
	max_threshold = amax(threshold*noise)
	    
	plot(self.wavelength, intensity,'b', label= 'time %2.1f ms' %self.time_stamp, linewidth=0.5)
	plot(self.wavelength,threshold*noise, 'g')
	


	xlabel('$\lambda$ [nm]',fontdict={'fontsize':12})
	
	if log_scale:
	    yscale('log', nonposy='clip')
	    axis([w_interval[0],w_interval[1],min_threshold/8, max(max_threshold,amax(intensity))])
	else:
	    axis([w_interval[0],w_interval[1],amin(-noise) , max(max_threshold,amax(intensity))])
	    
	
	legend(loc='best')

	
    def getData(self,sensitivity_correction = False,noiseCorrected = True):
	"""
	prepare intensities for  following processing. 
	In the fist step nonlinearity is removed, than if noiseCorrected is allowed the readout patterns, dark noise and hot pixels are removed and finally if  sensitivity_correction is allowed the different sentivity of whole device for different wavelengths is corrected.  
	
	"""
    
	#plot(self.intensity)
	#plot( self.Shot.readOutPatterns)
	#plot(self.Shot.darkNoise)
	#show()
	intensity = copy(self.intensity)

	if noiseCorrected:
	    
	    intensity -= self.Shot.readOutPatterns+self.Shot.darkNoise
	    intensity[int_(self.Shot.Spectrometer.hotPixels)] = median(intensity)
	
	
	correction = zeros(size(intensity))
	for i,a in enumerate(self.Shot.Spectrometer.linearityPolynomeCorrection[::-1]):
	    correction*= intensity
	    correction += a
	intensity = intensity/correction


	if sensitivity_correction:
	    intensity = self.Shot.Spectrometer.convertCountToPhotons(self.wavelength, intensity)

	return intensity

    def processSpectrum(self): 
	
	#if not self.plasma():
	    #return list()
	
	if len(self.lineFitList) != 0:
	    return self.lineFitList	
	
	corr_intensity = self.getData(False,True)
	#print 'identifyLines'#, self.Shot.readoutNoiseRMS
	self.lineFitList = self.Shot.Spectrometer.identifyLines(self.wavelength,corr_intensity,self.Shot.readoutNoiseRMS,  False ) # return table of lines, wavelength, +other parametres
	
	return self.lineFitList
  


     
 
def main(): 



    #print 'spektrometr'
    #x = linspace(0,1,500)
    #p = empty(5)
    #p[0] = 1
    #p[1] = 0.01
    #p[2] = 0.5
    #p[3] = 0.2
    #p[4] = 0.15
    
    #line = UniversalLineShape( p,x)
    #print lineCharacteristics(p)
    #print sum(line/500)
    
    #u = sum(x*line)/sum(line)
    #print sqrt(sum((x-u)**2*line)/sum(line))*2*sqrt(2*log(2))
    #print sum(x*line)/sum(line)
    #s = sqrt(sum((x-u)**2*line)/sum(line))
    #print (sum((x-u)**3*line)/sum(line))/s**3
    
    #plot(x,line)
    #show()
    #from SpectrometerControl import *
    #COMPASS = Tokamak('COMPASS')

    #COMPASS.saveConfig()
    ##exit()
    ##Spectrometer('gfdg',GOLEM)
    ##HR4000 = OceanOptics_Spec('HR4C1947',GOLEM)
    ##HR2000 = OceanOptics_Spec('HR+C1886',GOLEM)

    #GOLEM = Tokamak('GOLEM')

    #HR2000 = OceanOptics_Spec('HR+C1886',GOLEM)

    #HR2000 = OceanOptics_Spec('HR+C1103',COMPASS)
    

    #shot = HR2000.loadData('./data/', 'Java_Data_file', 'spectra') 
    #shot.saveData(noiseCorrected = True)
    
    
    #shot.plot_all( None, log_scale = True, sensitivity_correction = True)
    ##for spec in shot.spectra:
	##shot.spectra[0].intensity+= spec.intensity
    ##shot.readoutNoiseRMS*= sqrt(shot.n_spectra)
    ##shot.n_spectra = 1
    
    
    ##shot.tableOfResultes(details = True, sort = False, only_identified = True)
    ##shot.plotLinesTimeEvolution(file_type = None, only_identified = True, n_lines = 50)
    #shot.positionHistorogram()

    
    #exit()
    
   
    
    #shot = HR2000.loadData('./data_VW/','simple','spectra')  
    #shot.tableOfResultes(details = False, sort = True, only_identified = True)
    #shot.plot_all(None)
    #shot.plotLinesTimeEvolution('png',only_identified = True)
    #shot.spectra[2].plot(sensitivity_correction = True, log_scale = False, noiseCorrected = True, w_interval = None)
    #show()
    #hist(shot.spectra[1].intensity,70,  range = (-40, 100))
    #show()
    #exit()
    #shot.getData()
    #shot.plot( None, log_scale = False, sensitivity_correction = True )
    #plot(shot.darkNoise)
    #show()
    #shot = HR2000.loadData('./data/','SpectraSuite_HighSpeed','shotd') 
    
    
    
    
    
    
    
    
    
    #GOLEM = Tokamak('GOLEM')

    #HR2000 = OceanOptics_Spec('HR+C1886',GOLEM)  
    #shot = HR2000.loadData('./','Python_Data_file','spectra.txt_Data_')    
    #shot.plasmaShot()
    
    
    #shot.plot_all(None,  log_scale = False, sensitivity_correction = True, plasma = True)
    ##shot.tableOfResultes(details = True, sort = True, only_identified = False)
    #shot.plotLinesTimeEvolution('png',only_identified = True)
    ##show()
    #colour = ( 'g','r','c', 'm', 'y', 'k' )
    #shot.spectra[1].plot(  log_scale = False, sensitivity_correction = True)
    #for i, element in enumerate(GOLEM.elements):
	#for j,line in enumerate(element[0]):
	    #axvline(x = line , c = colour[(i%len(colour))])
    #show()
    
    
    #for i, element in enumerate(GOLEM.elements):
	#shot.spectra[1].plot(  log_scale = False, sensitivity_correction = True)

	#for j,line in enumerate(element[0]):
	    #axvline(x = line , c = colour[(i%len(colour))])
	#title(element[1])
	#show()
	
	

    #exit()
    
    
    
    
    
    
    
    
    
    
    
    
    #COMPASS = Tokamak('COMPASS')

    #HR2000 = OceanOptics_Spec('HR+C0732',COMPASS)  
    #shot = HR2000.loadData('./data_nova/02.12.2011/','SpectraSuite_xml','shot000') 
    
    #for i in range(1,shot.n_spectra):
	#shot.spectra[0].intensity += shot.spectra[i].intensity 
    #shot.spectra[0].intensity/= shot.n_spectra
    #shot.readoutNoiseRMS/= sqrt(shot.n_spectra)
    #shot2 = HR2000.loadData('./data_nova/02.12.2011/','SpectraSuite_xml','shot001') 

    #shot2.spectra[0].intensity = shot.spectra[0].intensity
    
    #shot2.plot_all( None, log_scale = False, sensitivity_correction = True)
    #shot2.tableOfResultes(details = True, sort = True, only_identified = False)
    #shot.saveData(noiseCorrected = True)

    
    #exit()
    
    #shot.spectra[0].intensity = random.randn(HR2000.resolution)*shot.readoutNoiseRMS
    #plot(shot.spectra[0].intensity)
    #show()
    #shot.tableOfResultes(details = True, sort = True, only_identified = False)
    #shot.plot( None, log_scale = True, sensitivity_correction = True)

    #exit()
    
    #shot.plot_all( None, log_scale = True, sensitivity_correction = True)
    #exit()
    #s = spectrometer.fittingFunction(array((1,0.1,0.5)),x)
    #Shot( 1200,  HR4000)
    #print strptime('Thu Feb 04 16:27:51 CET 2010', "%a %b %d %H:%M:%S %Z %Y")
    #shot = HR2000.loadData('./data_nova/4.02.10/', 'SpectraSuite_txt', 'HSA100610') 
    #figure()
    #shot = HR2000.loadData('./data_nova/4.02.10/', 'SpectraSuite_txt', 'HSA10069') 
    #print strptime('Thu Feb 04 16:27:51 CET 2010', "%a %b %d %H:%M:%S %Z %Y")


    
    #exit()
    #shot = HR2000.loadData('./data_nova/10.10.2011/','SpectraSuite_xml','shot011')  

    #shot.shot_num = hash('HSA1011a0')
    #print shot.shot_num
    #spec = shot.spectra[0]
    #spec.processSpectrum()
    #shot.processShot()

    #spec.plot()
    #shot.saveData(noiseCorrected = True)
    #shot.plot(None, log_scale = False, sensitivity_correction = True)
    #shot.plot_all(None, log_scale = False, sensitivity_correction = False)

    
    #shot.tableOfResultes(details = True, sort = True, only_identified = True)
    
    #sh993a
    #shot.processShot()
    #shot.tableOfResultes(details = False, sort = True, only_identified = True)
    #spectrometer.convertCountToPhotons(arange(200,1200),1,True)
    #plot(x,s)
    #show()
    #exit()
    
    
       
    
    #path = './data/'
    #COMPASS = Tokamak('COMPASS')
    #HR2000 = OceanOptics_Spec('HR+C1103',COMPASS)
    
    #HR2000.integTime = 1.8#[ms]
    #print COMPASS.actualShotNumber()
    #HR2000.start()

    
    #try:
	#dir_list = os.listdir(path)
	#for i in dir_list:
	    #os.remove(path+i)
    #except OSError:
	#os.mkdir(path)
	#print 'can not remove old files'
    

    #while True:
	#print 'waiting on software trigger'
	#COMPASS.waitOnSoftwareTrigger()			#TODO otestovat
	#print 'software trigger accepted'
	
	#print 'waiting on hardware trigger'	
	##COMPASS.waitOnHardwareTrigger(path+'spectra.txt')
	#print 'hardware trigger accepted'
	
	#shot = HR2000.loadData(path,'Java_Data_file','spectra')
	#shot.saveData(filename = None, noiseCorrected = True)

	##shot.plot_all('png', log_scale = True, sensitivity_correction = True)
	##shot.plot('png', log_scale = False, sensitivity_correction = True, w_interval = [200, 900])
	#shot.tableOfResultes(details = False, sort = True, only_identified = False)
	##shot.plotLinesTimeEvolution('png',only_identified = True)
	#COMPASS.storeData(path) 
	
    
	#exit()
    
    
    path = './data/'
    GOLEM = Tokamak('GOLEM')
    HR2000 = OceanOptics_Spec('HR+C1886',GOLEM)
    
    HR2000.integTime = 1.62#[ms]
    HR2000.start()

    

    

    while True:
	try: 
	    dir_list = os.listdir(path)
	    for i in dir_list:
		os.remove(path+i)
	except OSError:
	    os.mkdir(path)
	    print 'can not remove old files'
	
	print 'waiting on software trigger'
	GOLEM.waitOnSoftwareTrigger()			#TODO otestovat
	print 'software trigger accepted'
	
	print 'waiting on hardware trigger'	
	GOLEM.waitOnHardwareTrigger('ready')
	print 'hardware trigger accepted'
	    
	shot = HR2000.loadData(path,'Java_Data_file','spectra')

	shot.saveData(filename = None, noiseCorrected = True)

	shot.plot_all('png', log_scale = True, sensitivity_correction = True)
	shot.plot('png', log_scale = False, sensitivity_correction = True, w_interval = [200, 900])
	shot.tableOfResultes(details = False, sort = True, only_identified = True)
	shot.tableOfResultes(details = True,sufix = 'full', sort = True, only_identified = False)

	shot.plotLinesTimeEvolution('png',only_identified = True)
	GOLEM.storeData(path) 
	print 'data stored '+asctime(localtime())
    print 'success !!!!!'
    
    
    exit()
    
    
    
    
    #exit()
    
    
    #path = './data_nova/'   
    
    #GOLEM = Tokamak('GOLEM')
    
    ##HR4000 = OceanOptics_Spec('HR4C1947',GOLEM)
    ##HR2000 = OceanOptics_Spec('HR+C1886',GOLEM)
    ##HR2000 = OceanOptics_Spec('HR+C0732',GOLEM)    
    ##HR2000 = OceanOptics_Spec('HR+C1834',GOLEM)
    ##HR2000 = OceanOptics_Spec('HR+C1103',GOLEM)
    

    
    
    #working_dirs = os.listdir(path) 
    #database = list()
    #for directory in working_dirs:
	#shots = os.listdir(path+directory) 
	#for data_file in shots:
	    #if data_file[-8:] == 'SpectraSuite_xml':
		#data_file = data_file[:-9]
		#name = path+directory+'/','SpectraSuite_xml',data_file
		#print name
		#try:
		    #shot = HR2000.loadData(path+directory+'/','SpectraSuite_xml',data_file)
		    #shot.tableOfResultes(details = True, sort = True, only_identified = True)
		    #shot.plot_all(None, log_scale = False,sensitivity_correction = False)
		    #shot.saveData(filename = None, noiseCorrected = True)
		#except Exception as inst:
		    #print 'problem', inst
		    

    
    
    
    
    
    

    path = './data_nova/'   
    
    COMPASS = Tokamak('COMPASS')

    HR2000 = OceanOptics_Spec('HR+C0732',COMPASS)  

    
    working_dirs = os.listdir(path) 
    database = list()
    integtime = list()
    def ReadFile(directory,data_file):	
	name = path+directory+'/','txt',data_file
	#print name

	#if data_file[-3:] == 'txt':
	    #data_file = data_file[:-4]
	    #print name
	    #return 
	if data_file[-8:] == 'ProcSpec':
	    data_file = data_file[:-9]
	    #name = path+directory+'/','ProcSpec',data_file
	    #print name
	    

	    try:
		#print locale.getdefaultlocale()

		shot = HR2000.loadData(path+directory+'/','SpectraSuite_xml',data_file)
		#shot.tableOfResultes(details = True, sort = True, only_identified = True)
	    #if shot.readoutNoiseRMS != 0:
		integtime.append(shot.integ_time)
		if shot.integ_time > 6:
		    return 
		shot.plot_all(None, log_scale = True,sensitivity_correction = True)
		#shot.spectra[0].plot()
		#clf()
		#show()
	    except Exception as inst:
		print 'problem', inst
		    
	
	
    for directory in working_dirs:
	shots = os.listdir(path+directory) 
	
	#if directory!= '2012':
	    #continue
	for data_file in shots:
	    if os.path.isdir(path+directory+'/'+data_file):
		shots_2 = os.listdir(path+directory+'/'+data_file) 
		for data_file_2 in shots_2:
		    #print path+directory+'/'+data_file+'/'+data_file_2
		    ReadFile(directory+'/'+data_file,data_file_2)
	    else:
		ReadFile(directory,data_file)

    
    
    
    #integtime = array(integtime)
    #hist(integtime,50, range = (0, 20))
    #show()
    
    
    
    
    
if __name__ == "__main__":
    main()
    	
    	
    	
# vrací to oto blbě výslednou tlouštku    	
#TODO otestovat zabranení vícenásobného spousteni
# dívat tu javu v nekonečné smyčce? 
