From 7ec5b2968017e1d25e9605522afd6f5fbeae3fa4 Mon Sep 17 00:00:00 2001
From: Olivier Sauter <olivier.sauter@epfl.ch>
Date: Wed, 26 Sep 2001 09:21:10 +0000
Subject: [PATCH] add files needed for sxr and ECE in TCV/private mainly

git-svn-id: https://spcsvn.epfl.ch/repos/TCV/gdat/trunk@1823 d63d8f72-b253-0410-a779-e742ad2e26cf
---
 TCV/VsxrTCVradius.m            |  24 ++
 TCV/ece_te.m                   |  67 +++++
 TCV/get_xtomo_data.m           | 497 +++++++++++++++++++++++++++++++++
 TCV/private/angular_fact_1.mat | Bin 0 -> 3502 bytes
 TCV/private/bad_channels.mat   | Bin 0 -> 28504 bytes
 TCV/private/ece_d.m            |  51 ++++
 TCV/private/tcv_vesc1.mat      | Bin 0 -> 3285 bytes
 TCV/xtomo_geometry.m           | 402 ++++++++++++++++++++++++++
 8 files changed, 1041 insertions(+)
 create mode 100644 TCV/VsxrTCVradius.m
 create mode 100755 TCV/ece_te.m
 create mode 100755 TCV/get_xtomo_data.m
 create mode 100755 TCV/private/angular_fact_1.mat
 create mode 100755 TCV/private/bad_channels.mat
 create mode 100755 TCV/private/ece_d.m
 create mode 100644 TCV/private/tcv_vesc1.mat
 create mode 100644 TCV/xtomo_geometry.m

diff --git a/TCV/VsxrTCVradius.m b/TCV/VsxrTCVradius.m
new file mode 100644
index 00000000..38f1ce9f
--- /dev/null
+++ b/TCV/VsxrTCVradius.m
@@ -0,0 +1,24 @@
+function  radius=VsxrTCVradius(ymag,xchord,ychord)
+
+% get intersection of magnetics axis with chords
+% input:
+% Vsawtooth.shot : shot number
+% Vsawtooth.expdata.zmag : positions of the magnetic a 
+%        xchord:  two x-coordinates.
+%        ychord:  two y-coordinates. 
+%                 For each line (2xnl), they specify start + end pointsxis
+
+% Output :
+% Vsxr.radius : intersectionof magnetics axis with chords
+
+%parameter: start and end of each of the camera number 2
+
+
+% calculation intersections
+xchord=xchord/100;
+ychord=ychord/100;
+for i=1:size(xchord,2)
+  a=(ychord(1,i)-ychord(2,i))/(xchord(1,i)-xchord(2,i));
+  b=ychord(2,i)-a*xchord(2,i);
+  radius(:,i)= (ymag-b)/(a+eps);
+end
diff --git a/TCV/ece_te.m b/TCV/ece_te.m
new file mode 100755
index 00000000..32f05aa0
--- /dev/null
+++ b/TCV/ece_te.m
@@ -0,0 +1,67 @@
+function [TE_ECE,TE_ECE_ERR,RHO,R,T,TE_THOM,TE_THOM_ERR,Fcentral,CAL,NE_THOM]=ece_te(shot,Tc,T,SS)
+
+%	[TE_ECE,TE_ECE_ERR,RHO,R,T,TE_THOM,TE_THOM_ERR,Fcentral,CAL,NE_THOM]=ece_te(shot,Tc,T,SS)
+%	[TE_ECE,TE_ECE_ERR,RHO,R,T,TE_THOM,TE_THOM_ERR,Fcentral,CAL,NE_THOM]=ece_te(19314,[0.1 0.29],[0 1.5],10);
+%	
+%	Program that compute TE_ECE and TE_THOM profile over R or RHO
+%	for all the times T with the error matrices TE_ECE_ERR and TE_THOM_ERR.
+%	The calibration matrice CAL comes from the ece_calib.m program.
+%	Fcentral is the central frequency of the ECE working channels.		
+%	
+%	shot 	= shot number
+%	Tc	= [a b] temporal boundary on which we want the calibration on 
+%		thomson temperature. If you have no idea, put Tc=[0.1 0.29]
+%	T	= either
+%		-temporal vector
+%		-if T=10 => takes all the ECE times
+%		-if T=[c d];=> c and d are the temporal boundary for the output
+%	SS >= 1 	=> sampling rate SS if T=10 or T=[c d]
+%
+%	B est calcule avec la routine BandBres_allt de O.Sauter
+%	
+%	Blanchard  25.11.2000
+
+
+%-----------------------------------------------------------------------
+%	Recherche de la configuration de l'ECE
+%-----------------------------------------------------------------------
+[RHOece,Rece,Zece,Tece,Fcentralrho]=ece_rho(shot,T);
+[CAL,Fcentral]=ece_calib(shot,Tc);
+
+%-----------------------------------------------------------------------
+%	Recherches des signaux voulus
+%-----------------------------------------------------------------------
+[ECE,TECE,Fcentral]=ece_raw_signals(shot,T,SS);
+%[RHOece,Rece,Zece,Tece, Fcentralrho]=ece_rho(shot,T);
+%Zece=Zece*ones(size(Rece));
+[TEthom,NEthom,TEerr,NEerr,Tthom,RHOthom]=thom_rho(shot,6);
+if exist('TEthom')==0|length(TEthom) ==0
+   disp('Le profil Thomson proffit n''existe pas. On prend le profil direct') %(Tthom,Nbre de pts selon rho)
+
+  [TEthom,NEthom,TEerr,NEerr,Tthom,RHOthom]=thom_rho(shot,3); %(Tthom,Nbre de pts selon rho)
+end
+[t2,i3,i4]=common_ece(Fcentralrho,Fcentral);
+Rece=Rece(i3,:);Zece=Zece(i3,:);RHOece=RHOece(i3,:);
+TE_ECE=ECE.*repmat(CAL(2,:),length(TECE),1);
+T=repmat(TECE,1,length(Fcentral));
+if length(Tece)==1
+  RHO=repmat(RHOece,1,length(TECE));RHO=RHO';	
+  R=repmat(Rece,1,length(TECE));R=R';
+  X=size(RHOthom);
+else
+  %RHO=interp2(repmat(Fcentral,1,length(Tece))',repmat(Tece,length(Fcentral),1)', ...
+  %	    RHOece',repmat(Fcentral',length(TECE),1)',T');
+   RHO=interp2(repmat(Fcentral,1,length(Tece))',repmat(Tece',length(Fcentral),1)', ...
+	    RHOece',repmat(Fcentral',length(TECE),1)',T');
+  RHO=RHO';		
+  R=interp2(repmat(Fcentral,1,length(Tece))',repmat(Tece',length(Fcentral),1)', ...
+	    Rece',repmat(Fcentral',length(TECE),1),T);
+end
+ 
+X=size(RHOthom);
+TE_THOM=griddata(repmat(Tthom,1,X(2)),RHOthom,TEthom,T,RHO);
+NE_THOM=griddata(repmat(Tthom,1,X(2)),RHOthom,NEthom,T,RHO);
+TE_THOM_ERR=griddata(repmat(Tthom,1,X(2)),RHOthom,TEerr,T,RHO);
+TE_ECE_ERR=CAL(4,:);T=T(:,1);
+
+
diff --git a/TCV/get_xtomo_data.m b/TCV/get_xtomo_data.m
new file mode 100755
index 00000000..3652dd33
--- /dev/null
+++ b/TCV/get_xtomo_data.m
@@ -0,0 +1,497 @@
+function [sig,t]=get_xtomo_data(shot,t1,t2,dt,fans,angfact,tag);
+
+%                                                                  -
+%[sig,t]=get_xtomo_data_m5(shot,t1,t2,dt,fans,angfact,tag);           
+%
+%	INPUT:
+%       shot:	   TCV shot
+%       t1:      start time
+%       t2:      stop time
+%       dt:	     timestep
+%       fans:    camera switch, e.g. [0 0 0 0 0 0 1 0 1 0];
+%       angfact: relative etendue, size:	[20 x 10]
+%                tag	'full' or 'brief', indicates if it is
+%                a slow or a fast Pentland acquisition
+%
+% OUTPUTS:
+%       sig:	    xtomo signals, size:	[sum(fans) x length(t)]
+%       t:       times
+%      
+%	ATTENTION: length(time) may be shorter than foreseen !!      
+%
+% This routine works on Matlab5.
+% Original routine for Matlab4 by Anton Mathias.  
+%
+% Last update: 25-08-1999
+%
+%-------------MAC:[FURNO.MATLAB5.XTOMO]----------------------------------
+
+
+
+%---- get  data and offsets ---------------------------------------------------
+
+	% set a flag if flattenign of noise 
+	% for channels with low signal
+
+
+	iflat=0;
+	minsiglevel=0.01;
+	satlevel=9.9;
+
+
+	st1=sprintf('%6.4f',t1);
+	st2=sprintf('%6.4f',t2);
+        sdt=num2str(dt);
+
+%	sdt=sprintf('%6.4f',dt); % modified from old version             
+
+      
+	if shot >= 13836 & shot <= 13848 
+		tstart=-0.02;tstop=-0.01;		%  this  one is to be used only for
+	else						%  shot=13836 to shot=13848
+        	tstart=-0.04;tstop=-0.01;
+        end
+	
+	ststart=num2str(tstart);
+	ststop=num2str(tstop);
+
+if nargin <=6 
+
+	shot=mdsopen('eltca1::tcv_shot',shot);
+	t=mdsdata(['dim_of(\base::xtomo:array_001[',st1,':',st2,':',sdt,',*])']);
+	S1=[];
+	S2=[];
+	S3=[];
+	S4=[];
+	S5=[];
+	S6=[];
+	S7=[];
+	S8=[];
+	S9=[];
+	S10=[];
+
+	if isempty(t)
+		disp('get_xtomo_data: sorry, nothing to be found for this shot ...')
+		return
+	else
+		disp('*-------------------------------------------*')
+		disp('| get_xtomo_data: getting data from MDS     |')
+		disp('*-------------------------------------------*')
+	end
+                                                                                               
+	if shot>6768
+
+		if fans(1)
+			S=mdsdata(['\base::xtomo:array_001[',st1,':',st2,':',sdt,',0:19]'])';
+       			offset=mdsdata(['\base::xtomo:array_001[',ststart,':',ststop,',0:19]']);
+        		S1=S-repmat(mean(offset),length(t),1)';                       
+		end                                
+
+		if fans(2)
+		        S=mdsdata(['\base::xtomo:array_002[',st1,':',st2,':',sdt,',0:19]'])';
+       			offset=mdsdata(['\base::xtomo:array_002[',ststart,':',ststop,',0:19]']);
+        		S2=S-repmat(mean(offset),length(t),1)';                       
+	      	end
+
+		if fans(3)
+			S=mdsdata(['\base::xtomo:array_003[',st1,':',st2,':',sdt,',0:19]'])';
+       			offset=mdsdata(['\base::xtomo:array_003[',ststart,':',ststop,',0:19]']);
+        		S3=S-repmat(mean(offset),length(t),1)';                       
+	      	end
+
+		if fans(4)
+			S=mdsdata(['\base::xtomo:array_004[',st1,':',st2,':',sdt,',0:19]'])';
+       			offset=mdsdata(['\base::xtomo:array_004[',ststart,':',ststop,',0:19]']);
+       	 		S4=S-repmat(mean(offset),length(t),1)';                       
+		end
+
+		if fans(5)
+			S=mdsdata(['\base::xtomo:array_005[',st1,':',st2,':',sdt,',0:19]'])';
+       			offset=mdsdata(['\base::xtomo:array_005[',ststart,':',ststop,',0:19]']);
+        		S5=S-repmat(mean(offset),length(t),1)';                       
+		end
+
+		if fans(6)
+			S=mdsdata(['\base::xtomo:array_006[',st1,':',st2,':',sdt,',0:19]'])';
+       			offset=mdsdata(['\base::xtomo:array_006[',ststart,':',ststop,',0:19]']);
+        		S6=S-repmat(mean(offset),length(t),1)';                       
+	      	end
+	
+		if fans(7)
+			S=mdsdata(['\base::xtomo:array_007[',st1,':',st2,':',sdt,',0:19]'])';
+       			offset=mdsdata(['\base::xtomo:array_007[',ststart,':',ststop,',0:19]']);
+        		S7=S-repmat(mean(offset),length(t),1)';                       
+	        end
+	
+		if fans(8)
+			S=mdsdata(['\base::xtomo:array_008[',st1,':',st2,':',sdt,',0:19]'])';
+		       	offset=mdsdata(['\base::xtomo:array_008[',ststart,':',ststop,',0:19]']);
+        		S8=S-repmat(mean(offset),length(t),1)';                       
+	        end
+
+			if fans(9)
+			S=mdsdata(['\base::xtomo:array_009[',st1,':',st2,':',sdt,',0:19]'])';
+       			offset=mdsdata(['\base::xtomo:array_009[',ststart,':',ststop,',0:19]']);
+        		S9=S-repmat(mean(offset),length(t),1)';                       
+	      	end
+
+		if fans(10)
+			S=mdsdata(['\base::xtomo:array_010[',st1,':',st2,':',sdt,',0:19]'])';
+       			offset=mdsdata(['\base::xtomo:array_010[',ststart,':',ststop,',0:19]']);
+       			S10=S-repmat(mean(offset),length(t),1)';                       
+		end
+
+
+		sig=-[S1;S2;S3;S4;S5;S6;S7;S8;S9;S10];   
+		[satrow,satcol]=find(abs(sig)>satlevel);
+
+		if ~isempty(satcol)
+			i_tlimit=min(satcol(:));
+			if i_tlimit>1
+				sig=sig(:,1:i_tlimit-1);
+				t=t(1:i_tlimit-1);
+				disp(['get_xtomo_data WARNING: some channels saturated',......
+				', t2 changed to ',sprintf('%6.4f',max(t))]);
+	       		else
+		      		sig=[];
+				t=[];
+				disp('get_xtomo_data WARNING: saturations, no data returned');
+				return
+	       		end
+		end
+
+
+		angfact=angfact(:,find(fans));
+		angfact=angfact(:);
+		gains=get_xtomo_gains(shot);
+		% earlier than shot xxx:
+		perm=[1:180,[180+[2,1,4,3,6,5,8,7,10,9,12,11,14,13,16,15,18,17,20,19]]];
+		gains=gains(perm);
+		iact=find(fans);
+		nact=length(iact);
+		idetec=[];
+		for k=1:nact
+			idetec=[idetec,(iact(k)-1)*20+1:iact(k)*20];
+		end
+		gains=gains(idetec);
+		gains=gains(:);
+		fac=ones(size(gains))./gains;
+
+
+		sig=-sig.*(repmat(fac,1,length(t)));
+
+		if iflat
+			for ii=1:length(t)
+				imini=find( sig(:,ii) < minsiglevel*max(sig(:,ii)));
+				sig(imini,ii)=zeros(size(imini));
+			end			
+		end
+
+        	t=t';
+       		sig=sig.*(repmat(angfact,1,length(t)));                                           
+                                                                
+
+	else
+
+		if fans(9)
+			S=mdsdata(['\base::xtomo:array_001[',st1,':',st2,':',sdt,',0:19]'])';
+       			offset=mdsdata(['\base::xtomo:array_001[',ststart,':',ststop,',0:19]']);
+               		V=S-repmat(mean(offset),length(t),1)';                       
+	        end
+		if fans(7)
+			S=mdsdata(['\base::xtomo:array_002[',st1,':',st2,':',sdt,',0:19]'])';
+       			offset=mdsdata(['\base::xtomo:array_002[',ststart,':',ststop,',0:19]']);
+        		H=S-repmat(mean(offset),length(t),1)';                      
+       		end
+
+		%---- read gains of the current shot ---------------------------------------
+	
+		gains=get_xtomo_gains(shot);
+		g1=gains(1:20)';
+		g2=gains(21:40)';
+
+       		%------ solid angle factors -----------------------
+
+		aom1=angfact(:,9);
+		aom2=angfact(:,7);
+
+
+		%-----  calculate correcting factors ---------------------------------------
+
+		fac1=(aom1./g1);
+		fac2=(aom2./g2);
+
+		%--- result: offset-, calibration- and angular-factor corrected data -------
+
+		V=-V.*(repmat(fac1,1,length(t))); % vertical camera:
+		H=-H.*(repmat(fac2,1,length(t))); % horizontal camera:
+
+
+		sig=[H;V];
+
+        	t=t';                                              
+                                                                              
+	end
+
+else
+
+        if (~strcmp(tag,'brief') & ~strcmp(tag,'full'))
+        	disp('Only full or brief accepted as tag')
+          	return	
+	end                                 
+	
+
+	trace_tree=['\atlas::t_rex3_' tag ':'];        
+	trace=[trace_tree 'signal_inp:p'];                                         
+	%shot=mdsopen('eltca1::tcv_shot',shot);             
+ shot=mdsopen(shot);   
+	tref=mdsdata(['_tref=dim_of(' trace,'001)']);  %trace :p001      
+        dtref=mdsdata('_dtref = _tref[1]-_tref[0]');    
+             
+                                                                      
+	S1=[];
+	S2=[];
+	S3=[];
+	S4=[];
+	S5=[];
+	S6=[];
+	S7=[];
+	S8=[];
+	S9=[];
+	S10=[];
+
+        t_1=mdsdata('_t = dim_of(\atlas::t_rex3_full:signal_inp:p001)');
+        dt1=mdsdata('_dt = _t[1]-_t[0]');                    
+	npts=floor(mdsdata(['_npts=(',ststop,'-(',ststart,'))/_dt']));
+
+      
+                        
+if isempty(t_1)
+		disp('get_xtomo_data: sorry, nothing to be found for this shot ...')
+		return       
+	else
+		disp('*--------------------------------------------*')
+		disp('| getting Pentland acquisition: getting      *')
+ if strcmp(tag,'brief')
+    disp('|         BRIEF data from MDS                *');
+	else
+		  disp('|         FULL   data from MDS               *');
+	end
+		disp('*--------------------------------------------*');
+	end
+
+
+                             
+	if ((dt - dtref) < 1e-7 | floor(dt/dtref) ==0)
+		step='1';
+	else 
+		step =num2str(dt/dtref);
+	end 
+
+	
+	ind1=find( abs(tref-t1) == min(abs(tref-t1)));   
+        npts1=floor(mdsdata(['_npts1=(',st2,'-(',st1,'))/_dtref']));
+        sind1=num2str(ind1);                                                              
+        sind2=num2str(npts1+ind1);                             
+        
+	t=mdsdata(['data(_tref)[',sind1,':',sind2,':', step, ']']); 
+                             
+        
+
+                                                                                      
+	if fans(1)
+		channel=['001';'002';'003';'004';'005';'006';'007';'008';'009';'010';...
+			'011';'012';'013';'014';'015';'016';'017';'018';'019';'020'];          
+		for i=1:20
+%			S(i,:)=mdsdata([trace, channel(i,:),'[',st1,':',st2,':',sdt,']'])';
+
+                       	
+			S(i,:)=mdsdata(['data(',trace,channel(i,:),')[',sind1,':',sind2,':' , step,']'])'; 
+       	   	     	offset=mdsdata(['data(\atlas::t_rex3_full:signal_inp:p',channel(i,:),')[0:',num2str(npts),']' ]);
+			S1(i,:)=S(i,:)-repmat(mean(offset),length(t),1)';                       
+		end                                
+         end
+
+        
+
+        if fans(2)
+		channel=['021';'022';'023';'024';'025';'026';'027';'028';'029';'030';...
+			 '031';'032';'033';'034';'035';'036';'037';'038';'039';'040'];
+		for i=1:20
+
+%			S(i,:)=mdsdata([trace, channel(i,:),'[',st1,':',st2,':',sdt,']'])';
+                                                
+			
+			S(i,:)=mdsdata(['data(',trace,channel(i,:),')[',sind1,':',sind2,':', step, ']'])'; 
+      	   	       	offset=mdsdata(['data(\atlas::t_rex3_full:signal_inp:p',channel(i,:),')[0:',num2str(npts),']' ]);
+			S2(i,:)=S(i,:)-repmat(mean(offset),length(t),1)';                       
+		end                                
+         end                   
+         if fans(3)
+		
+		channel=['041';'042';'043';'044';'045';'046';'047';'048';'049';'050';...
+			'051';'052';'053';'054';'055';'056';'057';'058';'059';'060'];
+		for i=1:20
+%			S(i,:)=mdsdata([trace, channel(i,:),'[',st1,':',st2,':',sdt,']'])';
+                                                
+			
+			S(i,:)=mdsdata(['data(',trace,channel(i,:),')[',sind1,':',sind2,':', step, ']'])';  
+        	       	offset=mdsdata(['data(\atlas::t_rex3_full:signal_inp:p',channel(i,:),')[0:',num2str(npts),']' ]);
+           		S3(i,:)=S(i,:)-repmat(mean(offset),length(t),1)';                       
+		end                                            
+         end
+
+        if fans(4)
+		channel=['061';'062';'063';'064';'065';'066';'067';'068';'069';'070';...
+			 '071';'072';'073';'074';'075';'076';'077';'078';'079';'080'];
+		for i=1:20
+%			S(i,:)=mdsdata([trace, channel(i,:),'[',st1,':',st2,':',sdt,']'])';
+
+                       
+			S(i,:)=mdsdata(['data(',trace,channel(i,:),')[',sind1,':',sind2,':', step, ']'])';  
+                      	offset=mdsdata(['data(\atlas::t_rex3_full:signal_inp:p',channel(i,:),')[0:',num2str(npts),']' ]);
+           		S4(i,:)=S(i,:)-repmat(mean(offset),length(t),1)';                       
+		end                                
+         end
+         if fans(5)
+		channel=['081';'082';'083';'084';'085';'086';'087';'088';'089';'090';...
+			 '091';'092';'093';'094';'095';'096';'097';'098';'099';'100'];
+		for i=1:20
+%			S(i,:)=mdsdata([trace, channel(i,:),'[',st1,':',st2,':',sdt,']'])';
+
+			
+			S(i,:)=mdsdata(['data(',trace,channel(i,:),')[',sind1,':',sind2,':',step,']'])';  
+                        offset=mdsdata(['data(\atlas::t_rex3_full:signal_inp:p',channel(i,:),')[0:',num2str(npts),']' ]);
+			S5(i,:)=S(i,:)-repmat(mean(offset),length(t),1)';                       
+		end                                
+         end
+
+        if fans(6)
+		channel=['101';'102';'103';'104';'105';'106';'107';'108';'109';'110';...
+		  	 '111';'112';'113';'114';'115';'116';'117';'118';'119';'120'];
+		for i=1:20
+%			S(i,:)=mdsdata([trace,channel(i,:),'[',st1,':',st2,':',sdt,':',step,']'])'
+   S(i,:)=mdsdata(['data(',trace,channel(i,:),')[',sind1,':',sind2,':', step, ']'])';  
+			offset=mdsdata(['data(\atlas::t_rex3_full:signal_inp:p',channel(i,:),')[0:',num2str(npts),']' ]);
+           		S6(i,:)=S(i,:)-repmat(mean(offset),length(t),1)';                       
+		end                                
+         end
+        if fans(7)
+		channel=['121';'122';'123';'124';'125';'126';'127';'128';'129';'130';...
+		 	 '131';'132';'133';'134';'135';'136';'137';'138';'139';'140'];
+		for i=1:20
+%			S(i,:)=mdsdata([trace, channel(i,:),'[',st1,':',st2,':',sdt,']'])';
+
+			
+			S(i,:)=mdsdata(['data(',trace,channel(i,:),')[',sind1,':',sind2,':',step,']'])';  
+			offset=mdsdata(['data(\atlas::t_rex3_full:signal_inp:p',channel(i,:),')[0:',num2str(npts),']' ]);
+           		S7(i,:)=S(i,:)-repmat(mean(offset),length(t),1)';                       
+		end                                
+         end
+        if fans(8)
+		channel=['141';'142';'143';'144';'145';'146';'147';'148';'149';'150';...
+		 	 '151';'152';'153';'154';'155';'156';'157';'158';'159';'160'];
+		for i=1:20
+%			S(i,:)=mdsdata([trace, channel(i,:),'[',st1,':',st2,':',sdt,']'])';
+
+			
+			S(i,:)=mdsdata(['data(',trace,channel(i,:),')[',sind1,':',sind2,':', step,']'])';  
+                        offset=mdsdata(['data(\atlas::t_rex3_full:signal_inp:p',channel(i,:),')[0:',num2str(npts),']' ]);
+           		S8(i,:)=S(i,:)-repmat(mean(offset),length(t),1)';                       
+		end                                
+         end
+         if fans(9)
+		channel=['161';'162';'163';'164';'165';'166';'167';'168';'169';'170';...
+		 	 '171';'172';'173';'174';'175';'176';'177';'178';'179';'180'];
+		for i=1:20
+%			S(i,:)=mdsdata([trace, channel(i,:),'[',st1,':',st2,':',sdt,']'])';
+
+			
+			S(i,:)=mdsdata(['data(',trace,channel(i,:),')[',sind1,':',sind2,':', step,']'])';  
+                        offset=mdsdata(['data(\atlas::t_rex3_full:signal_inp:p',channel(i,:),')[0:',num2str(npts),']' ]);
+			S9(i,:)=S(i,:)-repmat(mean(offset),length(t),1)';                       	
+	end                                
+         end
+         if fans(10)
+%		channel=['181';'182';'183';'184';'185';'186';'187';'188';'189';'190';...
+%		 	 '191';'192';'193';'194';'195';'196';'197';'198';'199';'200'];
+
+% this permutation takes into account some not well defined (not yet)
+% hardware corrections. Are the gains to be permutated ?
+
+
+		channel=['182';'181';'184';'183';'186';'185';'188';'187';'190';'189';...
+		 	 '192';'191';'194';'193';'196';'195';'198';'197';'200';'199'];
+		for i=1:20
+%			S(i,:)=mdsdata([trace, channel(i,:),'[',st1,':',st2,':',sdt,']'])';
+
+			
+			S(i,:)=mdsdata(['data(',trace,channel(i,:),')[',sind1,':',sind2,':', step,']'])';  
+                        offset=mdsdata(['data(\atlas::t_rex3_full:signal_inp:p',channel(i,:),')[0:',num2str(npts),']' ]);
+			S10(i,:)=S(i,:)-repmat(mean(offset),length(t),1)';                       
+		end                                
+         end
+
+	sig=[S1;S2;S3;S4;S5;S6;S7;S8;S9;S10];   
+
+
+	[satrow,satcol]=find(abs(sig)>satlevel);
+
+		if ~isempty(satcol)
+			i_tlimit=min(satcol(:));
+			if i_tlimit>1
+				sig=sig(:,1:i_tlimit-1);
+				t=t(1:i_tlimit-1);
+				disp(['get_xtomo_data WARNING: some channels saturated',......
+				', t2 changed to ',sprintf('%6.4f',max(t))]);
+	       		else
+		      		sig=[];
+				t=[];
+				disp('get_xtomo_data WARNING: saturations, no data returned');
+				return
+	       		end
+		end
+
+
+	angfact=angfact(:,find(fans));
+	angfact=angfact(:);
+	gains=get_xtomo_gains(shot);
+	% earlier than shot xxx:
+	perm=[1:180,[180+[2,1,4,3,6,5,8,7,10,9,12,11,14,13,16,15,18,17,20,19]]];
+	gains=gains(perm);
+	iact=find(fans);
+	nact=length(iact);
+	idetec=[];
+		for k=1:nact
+			idetec=[idetec,(iact(k)-1)*20+1:iact(k)*20];
+		end
+	gains=gains(idetec);
+	gains=gains(:);
+	fac=ones(size(gains))./gains;
+
+
+	sig=sig.*(repmat(fac,1,length(t)));
+
+		if iflat
+			for ii=1:length(t)
+				imini=find( sig(:,ii) < minsiglevel*max(sig(:,ii)));
+				sig(imini,ii)=zeros(size(imini));
+			end			
+		end
+
+        t=t';
+       	sig=sig.*(repmat(angfact,1,length(t)));                                           
+                                                                
+
+end
+
+mdsclose
+
+
+return
+
+
+
+
+
+
diff --git a/TCV/private/angular_fact_1.mat b/TCV/private/angular_fact_1.mat
new file mode 100755
index 0000000000000000000000000000000000000000..85fd9a52692c72745b87ad324e2601e6d90fa66d
GIT binary patch
literal 3502
zcmd5;i#wEQ8y{A(vKo?%ksW6t64tr<Movj(ZIgUf2}>nP&K<B)Ydc8Riky{FlEZ{$
z$svcKG)}d$7Hz&3m7-O)!+!Iu_uaj&{R_V5x@PWqp7%aI_x=0b@6;?|Fyt5vhNSRA
zmW1%*6CCiBk8c=b1FBva7ufD}z{IH<<)?0TSXG&C-_(g1^u+FDe3}(@zm^}=Jz<7}
z-6!Nlbxkm09k0W~d@Z_|UT({oHoz%2<;2fxbkR^^LUSNW3;UyH-#b*PBeB*HSH<78
zlbct0aWLP$u~Rl-C1wa$W(FHW136_`w^>vi7OIOI4Cg5-U_xWC?CEqC3RHKFC)pUF
zY(jeQmeV@aa}7M5KC;=*QWd$<3S!}@%IK9$A}CV+S0Kr9gOnVO3v70pCM`p4ujb#q
zl4bEhp|#}2{j{XdX+7GG_Cx!n<Ir*GJm|bwC?#6!Y~rneZZ)KWEGoVt?qc@3UVX_%
zxze$iQpc549XY6+NHkUP+2O#HvQTxbnkDtt!l8L`V_lrST#+rkdJUR>+~)G)@>&$V
z<}Br1Yl5ZiF>+30=9nF_x`Hchjc}8Aut70ys=|+9>#@p+_;JL;u_@0E2D+fKJ&9wB
zV-}^3hve2_;|AG;^TFmABj-3KS87bz8)BeLS!gRqkE*v8dYJEG_PMIzVIY1uNV+ky
znN2-kfu)Wl4>>gPCK|FxZhr)rnDQ(%P^FBCnubvmANrX{vWnfq#D?6C2y>nclJ1`m
zUWT3d35_AL@|a}fVtS>`35U$dz3q{y;JRz3*A`>GAo&4>H%lZ{YUP-s&8TToyNfYO
zUne&<#G{i#^F68tn5!7k!0XVZ?xlqiC7OoO6Y5m|Rk3it<bwm&9ID<c@n!GZ;(brp
z$n+rbSQK6r@HXetg)e?xfo(3N1M*n>^2ytiK{PH_Q1%Mcy;+#ttnJnEI~z$iSYpal
z-8d8uRH<_fEKwrN(57&ri)Ep=IYk!?C>*TCVRc1=hFoJj`83T=Jl+g{)p+eb*lmdh
z!?}u)4T$4Rf@510j~`ZN8`z_+BSD-KN(7O6%2V}ThMP}Yj(GaXP;r?Q|CmS^&tNdA
zJec^tziPI?UKYWG;71O7UCB*XP<&=n^`!f*M)69E>boBLSt@x-lo(QZ8dLQ&r{=;s
zim$fVu!Hp78T)I<JaE9HI;3BA6fY1fzLZ`q#<imGWrki$36G5_d>EojHJLB^fAind
zk5p2V;*mN|X!b|Fdar^6En_hp%I`|5uWTw`7LuxzDABxBpm55<1%EV;X<pMjr+H7$
z0eUXbc%bJ7JxA!dLQDFb)}!rcKeS&u4jq@ygU(A4=Z~rmjSE;9pGGDmi{d4V>bD}M
z2b0dTsrju;)q#Vg`b-a1s=w;g+||O9K?G+yDEVN&^@Bov3g?Cx`Wu<Y##qZId^W?k
ztAd`aKWT*-EHbZ=;-xKi%n`r#xaZycQMI9(f4y543$f-Cvd=%jhmj7O8oavyIDP@H
zB~Lr459EQ0gu#dV0|nrEqi#n2>q}tMU9j&P-(px4dA#Yd+jXd~8xHEsFNYBBwKVny
z9$4zSD-FJ>2DW#>rn&f9Fr3^TtkBc|FB>1JONBK-zVH5BQY9^LZ?}nF+Gnk>YTFWP
z59xMDm26omvZDhI_{6Taj^_iz`ht4HEk68`5H$0(4j)1|qhbMOJg}GA=2^M41~MG@
zMV8t15N(^UJL=d3wMUy14Y{{r@pYXle~S)aN{R<%=5|53*QwTBKX*e8*J9I;904S5
z=FYm_6u_ccqb6tef^-u~a{iy23wq}`3wCxtFZfB$UGTftVj+$^myWxT$LHA#d1Wl+
z!w0*8tg;h)@H+FV(7B5b*H(MpPe0oM7FLr#o;ubB{(ZOp+}_a&>}5{*J>xB4p<Fk>
zI?)7v%U=y06f}T{R`<`J_SHg=hs^Vf!&R{8rxywm{7OixSRSujTn436E>5*&*Fj4}
zrhWX~6$l;J!M^+95?I>p@oSzh03);SLuWd3p?1B+c5Ox`ydC5J0TJn-rek@>${`<O
z+|T3~BwvJM4<|nvl`V#!3|<+uXO)6d{C_7ZpI5-sxGDM6S5;sX_Utoby*j8fnwfuR
z)d(4ZAyW+Z7HBN`YD2qkJD56b%8EAZf)urx@BoE-pd|8id0A%<jM_~+ynDSD1Xl_w
zW!LmUQtUwe@{N7apL<0yFs2utz2BX?C$$Hxn>@|Jg9NZr{!P~4sXM^i!K)iq&IH9}
z$2<>3WJAh~>^Et=e3-K+X}FtJ2#rHsD_vI=f!bu5k?FaskRG0+u|~5L22bBrROu>*
zkvzq-VX{2nZEE+}7FrGc?m^y95^7=NjqgnQ4>dsb_68GUohFD=ad_c5+5)E9VWJ;p
zw!*E|VTES1?XZ~@DsDX00n57!Q$1ulK|_D!d1q7?>{ou#saMklb4fX#izK=r=*#u7
zw%a-(_3@u!-zD=wcIuQ~*rpDsPK?|hk<kW6?j=2wa=Hz7$3Im&@O3j>eSD)SrL_?j
zJq}=GzOIK0PlqN$Eovcl%eTj3$E#sTR7Jskga=0?*IEnO%Rw^G+wa`|Qg~<J2I(eO
zAx~S(d5J|4?7C*i)%PtF)>R<r{v{vUWWG&0ppy+cW#(Tv9nXN;VI~f7??7VM`nVU5
z1aMAsTEIQq16Lkw&Gz5h2OO5dMw_#JaHd|jpf<S=TC$>5TMqX^W%IG~VbMKs*dyq5
z%EUcTQC^ib|ELSPjUyjq541y+tE=;E?-mFsI(2?&ZzC-B4d>}=*MYvOue9J~6}asf
z92xte0v_#I_uG;AQaH1<aQ;wRF~o^$oY2+32+q7)en%wpp?=Mly`N;B2eo9iLYb#o
z;2G;LKkJtV(JAg4TecR0{YK5>UriUmuW^s;&Hb-}q|H6&twE(=oO5SO{LylF*q406
z{|yhODw}GgcT~e~)1K7=tvYaz=XmGM-GYj3O9r+ZG{M}}-uv9=E#Q0QyQyo&t?<){
z7L#h-c8Jt|-*j5916l)pc1C*eVOm#o=8h2`GCk}AEoM5Pd1TZoTj&3MzO+O0{qS_b
zs}9&(?)DnB_|PEM`p8GE19tw96d_X54jac0_IQ5W21|04VESz<a7EqTJhTzkL;vNn
zoR6Ep)+WMle0c-t6~p^@el66$9ukORHH3U)nP|9`2b?`^^WJ*pP+j1urF^Lrw)*z{
zVi+jIq54Ce__-o*IOHuQl28cx+xp}_vd@Fb%T+g3zt4j2>o#Zi-}*uLY7r5>PT1uC
zi*O0IjPT&yVF$xPcLxVBjvxcc@gZcy5g9}6h39{*!^mhPKXxMdW-y}u(QhPsp&b$#
w{f}%Uds=TFlHEV{v_35(skpoGFSWnp`2TB9=RxO1#}l5%5lcFr@cgg!KdxQ6Gynhq

literal 0
HcmV?d00001

diff --git a/TCV/private/bad_channels.mat b/TCV/private/bad_channels.mat
new file mode 100755
index 0000000000000000000000000000000000000000..e950698cfd82fae24b1fcd0aeaf8eea96ef638b1
GIT binary patch
literal 28504
zcmeI2X>3$g6vtoLK@izQ1$Ds0ParK=G>uz9k<PmV3I!JsEm8z6B~}m>`2ZsByFz5M
z2+=PtAB}Nea7p}NOh5$zKe$8%H5yk`{2kh3u8;fXy-rI@;ojuOxyw1{-uFLmrf-JH
zl~<bbF`3fSiZi+L%Hgx;)>IG6TvSuFVD|jFd1aZ3F}W)zO&*pRTUTARpn6tjer;K%
za$$8QyRarxQk*FnT~<1(Y}A>Vk;TO&nWna)AXgCt!O=lb`ysXqbZfMjpmW2xXePFQ
zP!trE%k2L7+WLL8hT*tvK`?GR?0}s?Fn$Ycg^%GA_!L+hFKfGkprR|VPepe)3VK2>
z=ncm}ALtAHpg#<Nli*|+0z=^xI2DFN4ko}vxELnEWT=EGFcmI`X)qmTz*SHM)i4{b
zfjKZ2u7&HM2Cj!XxB(WxLbwqY!A)>8+yaZ?R#*Z{VHvE1Rq!~hhBdGj*1>vs0_tG{
zJP8|N6TAQ~!>jNbybf=`oA4HFhIimycn{tWf@~3VfR4}!s9lzNX3;u})>*X9_5k$G
zqG=XQv&X_fI35N88fVcsJ2VKw6M(sJFq{ZOK<0-5nujBx7*2;$7!7B@S#UO-1Lr|G
zjDfLmK3o9f;6fM=6_AAxsAEVSL+Th(x9}3U6fOhm7fyrea0Sc&Y8Y0*Ot>0m0rd>2
zXGlH6IY1r5c~A@UfqI72Go+p&wG63cxCp3YxEPkeQdkDJ!E(4AR=^!_C)^JYz=QA*
zY=WoYX?O;ng|}fdd<0uyD{O=9@EPoYo$xt)0lQ%jd=EdskMI-x48Opy@EiONf52Y&
zD+mnrG(Djg91VS-FZ2WIY^b9d1cM<1)X<CoYG+DdB#Z*;Xi9-P8tQ1MpE(z(pP^=k
zni*<lsF|T=hME~_W(?$D0!)NSFa@Rp^)uAZP(MTc%ru~WhWZ)mXQ-K(1=P$?GsEl7
zP%l#h^Pm>$fLa;qW2lE&4%ER=2XhzP4fnvka34GjkHDkw7_5X<uo~9DT383`p&mBC
zMtBaMhZo@`cm=4R`2eV!*$&jqP&2a&sGIo`zJjme8~7H!gWa$PsG*^DhT0iwXQ-W_
zcIHpm3sS$o;UCx+1jrcVil75@gwD_fx<WVT4i$j*xiBbNQ6F^J$)3D3d%Onz9`N^o
zzX#g#9^idD=-o&!tFLp-v69cuspGU(&L*fo&kFm`8kbYg8eNw-drtkLR`(IxsS)vU
z#qu}qKixEs_{qHBic!Bq&FVd%dh7XmI<8UI;w*pfgXR#IH{ZNyZR>6KB~I5}mgkiJ
z#pc^b*P?lPkVDXQmyO26jU_g^#xvFaIkmFx+9MjPKK5%Un%gzjIM-UzoW`~%+mvfc
z@^!6?7WTYL8=jqm0gcGI_LD51nSvyAKlRgbS}SK0Y=1eIjP-ocwKVG*8?|;Wx#xvw
zK8>-@V9&|CXf9T-`(Ri;y9P>J{(<-;Avr{2@l(AuMy^BqX|1hcA^o)b*nYZ?&dFH&
zVs+hRwN|R7%i6hI`Sdjw^|QyWy<Ky0`NdbyDt$Y*-M4e<eIa#~XA#BrarHIMwU#ue
zu}RiGb4&H(_H@-$=Vb4sYf@jaNxrW2Xq<h$`raYWeL64x+E20?pU<i;>SJGz_H`fe
z(K)sHOsZq?5q}wLoUTRZ>~+npb8*(rrO$NYYv<NjJ&TUT%A6pzaK+fZuD-^nv8t`u
zIDeKBo18;*mZpDu8b%G*q5ZVhbzM!FOS_Nlr*+dg8QV?QT{cav>|D~1=F!())XyGA
zW8%hf`Ndz)D*d=-=T+UMzIsMGhwW?kHO`((+t=768+ZQZ%*VYB@|{KYOS(2YKFQa$
zF7fufeT|yio|EfH+cZ9(RbBG==zij({%ZAoLC4}FO~%@EEjnkfYi^y3vvw}MZ^hTn
zt+9F*9gCGYK{R*8FrBw)d>U(?HO@cpX4azUEX}@V8B5;w(DeHCtkHApy30zfqW+pE
z&PMa-Yb?%R=A$ujW0$X<QEk5S+IjUnQeQozSo?3dtFLjcwWK+XZ3}B(yWLMZi|y;`
z+wrdTq&bPVee8U$zWBCe)4kK?RDZSlU8{~YkBqhHT6E4{m%hEOYf)?G()%`Ttes2V
zBh&0>$27HXpS7v6QGa>=m8LlxMn}fKEm5a~t(%@vZBr{FuOs{YQ_mY`HGkZCzIlyp
zkM_9sc<5F*vwk<w3R@kpg+6Z{u#jRx?$h?v+S_^CleUM#S>^vC?f1?@;Y|Jv?MDIq
z|BtpfvEH|BvA5sTE6|b(@ViytyaEL)5dFTg;6y&qD{urWAfGFJ^9mHIfc+h9p*eh<
zSKvrjKt5Oc<`pPZf#!TiTj-uX&MR=ZDWLxkd$?Vfzn<1nz<)>EIy(DZ4sQkIbER)y
zfkG9K-;wy=k`$_*kMjyNC?LOC^UW*Z74Qmp1=>yl|GUJt)86mw74Qmp1-t@Y0k42p
sz$@St@CtYZyaHYUuYgy;E8rFI3U~#)0$u^HfLFjP;1%!+9HI*R3ydd$&j0`b

literal 0
HcmV?d00001

diff --git a/TCV/private/ece_d.m b/TCV/private/ece_d.m
new file mode 100755
index 00000000..17752b9d
--- /dev/null
+++ b/TCV/private/ece_d.m
@@ -0,0 +1,51 @@
+% 	Fonction matlab  donnant la relation entre la frequence centrale des 
+%	cannaux du radiometre avec le no du cannal
+D=[  78.8500   12.0000
+   79.6000   13.0000
+   80.3500   11.0000
+   81.1000   14.0000
+   81.8500   10.0000
+   82.6000   15.0000
+   83.3500    9.0000
+   84.1000   16.0000
+   84.8500    8.0000
+   85.6000   17.0000
+   86.3500    7.0000
+   87.1000   18.0000
+   87.8500    6.0000
+   88.6000   19.0000
+   89.3500    5.0000
+   90.1000   20.0000
+   90.8500    4.0000
+   91.6000   21.0000
+   92.3500    3.0000
+   93.1000   22.0000
+   93.8500    2.0000
+   94.6000   23.0000
+   95.3500    1.0000
+   96.1000   24.0000
+   96.8500   12.0000
+   97.6000   13.0000
+   98.3500   11.0000
+   99.1000   14.0000
+   99.8500   10.0000
+  100.6000   15.0000
+  101.3500    9.0000
+  102.1000   16.0000
+  102.8500    8.0000
+  103.6000   17.0000
+  104.3500    7.0000
+  105.1000   18.0000
+  105.8500    6.0000
+  106.6000   19.0000
+  107.3500    5.0000
+  108.1000   20.0000
+  108.8500    4.0000
+  109.6000   21.0000
+  110.3500    3.0000
+  111.1000   22.0000
+  111.8500    2.0000
+  112.6000   23.0000
+  113.3500    1.0000
+  114.1000   24.0000            ];
+  
diff --git a/TCV/private/tcv_vesc1.mat b/TCV/private/tcv_vesc1.mat
new file mode 100644
index 0000000000000000000000000000000000000000..becc16ea4e9dc15a3afedbacd112ef4f681b7bcb
GIT binary patch
literal 3285
zcmeH|e@v8h9LEo)SdDDW<cF)+NUelYnv}=<@exWnLAW2n@#9WT1vlD+A&Stz87^fq
z!4Q9-e_R=2a3GkJ;PK<bPCzUz9Q8*`wD=<tPP(wDK<a?cn>Aa1vC%)?_Uv=dbI<d9
z-|yG^z56N@3PqChaCIKzA2;W*^SzAZl!-Hi!bp=}W`1^9(%r9C8t7E;Pk|YI@l<3O
z3;A<VJRQoIuid*fj&hIqn_SJYG(R=*#Xmleq2r&qr&WHTCy(E1Z3&m6>4ehgF|RI)
zymDqIx7BT;lU!xZs@IYCtmwI~UDi^``Z3d;%NkmcUv}}RRZUm=7RDR*MbfqWaIRKG
zi0{X9GXH;b04Hz*M{osaa7PZv1vw!%<cM66Gjc}{(2M`8myfPHSU6m9tnaTiLlz}1
z8}uYqU0}gO$4gXnEWvHtfFHFt6s$-b3ZT~WcKh+41IhML_$M89C2^4@X7t`FN(i2P
z(DMP&W^3&aHNOOtae3m(4UV_TP&MyF@8CL$Uu<(^-w&lY^UiP&`v!{j@;z7GA5JkJ
zcDktUMUZXowwt#NYKpi>)({sBb-9<H$vmr}xb4p#s#0lb*d=twk9O^3pGi0Ed8@y)
z$=^>$1D*p1cBE*@vi4+mW4ea?KA_jrcB!XyK__&>2YkUNe1ikHfD^cZBe;SyxFZLd
z&2pXPK9=WL-ox@<mUURxWzoT+i$y1kZWbRbzF2&+<iWy$g$oNO7H&MzAJR>GO0Fz9
zv*gaw1D0M)>4Hw^h7b6HPxuB0Z~-T914nQLXK+Uj$YtF3KX;lr&U5EZle&A?_%8~D
zNw%07oLAj9wKeZClF2eF-{vsT{$#J!zONXlxFM*0-{yEKd&a1~T^vW>>@Cf6|1y?(
z3cHt@zm1`WCv?SOH}v%Wn$gkWwrEPIOW#u397QVw`dih`t;V}z>yn2xIvV!<zM#BZ
zOI@cmLs4fmbm6%Uz3XW;t*S4h?E%iMW-M#LvdWD#;h!GQ|Kb2n;0BK13eMn;9FPlg
zLT<<rxguxejvk;F=*fh>8SLky>)pzH-lE*CAjbvM!SAZQsJ5VDabD3v>S(|3y6Cz;
z-N?Jyb7%h>6npkT_|lKxq}~<}Yfe)TrSz<L^{ejHl$}!*b~t-2wd6PDF3%64Ij6p^
z?=E?Vvi)=J=6<@KT7n9EmRwSio5~#B*&0S+`kKR;T^s4BCdq40cO*5A^ftbGO--rm
zuIoFTfB8>F-u?CE7qnz8=`1|ZrX_ur>y~q=vc>T9=To60KBiI~op+QBm`AkqbW57|
zTM8}pEBkVKe^=9yq?Ri_y^&Kop&LHn3qIi+9KZ#fzzrP16`a8xIml$W&T=2ib1d&+
zc`wU4EbFr9V9~{*lSMa+4;Ei6d9e6q;lRR$g%b-mp6CzDRx^z$xw7QUk~>QeSbD+I
elPR6h4Il6YpYRP1-~vwI29Dqg&ftz5kjvj<o-(!o

literal 0
HcmV?d00001

diff --git a/TCV/xtomo_geometry.m b/TCV/xtomo_geometry.m
new file mode 100644
index 00000000..a809b9e6
--- /dev/null
+++ b/TCV/xtomo_geometry.m
@@ -0,0 +1,402 @@
+function [fans,vangle,xchord,ychord,aomega,angfact]=xtomo_geometry(i_detec,fans)
+
+% ----[anton.public]
+%
+% function 
+%
+% [fans,vangle,xchord,ychord,aomega,angfact]=xtomo_geometry(i_detec,fans);
+% 
+%	inputs: 
+%
+%  i_detec:	=2: Xtomo prototype cameras (shot# < 6768)
+%		=1: Xtomo 9-cameras	    (shot# > 682x)
+%
+%	outputs:
+%
+%  fans:    camera switch, 1=on,0=off (1x10)
+%  vangle:  angle between detect. surface normal and pos. x-axis (1x10)
+%  xchord:  two x-coordinates (2xnl) and
+%  ychord:  two y-coord. for each line (2xnl), they specify start + end points
+%  aomega:  etendue in mm^2 x steradians
+%  angfact: angular factors, inverse of relative etendue (throughput) (20x10)
+%
+%	uses:
+%		AOMEGA=etendue_n2(b1x,b1y,b1z,b2x,b2y,b2z,z01,z02,X0,cw);
+%		angular_fact_*.mat , '*'=i_detec
+%  
+%
+%---------------- M.Anton 14/3/95 -------------------------------------------
+
+disp('*----------------------------*')
+disp('|   this is xtomo_geometry   |')
+disp('*----------------------------*')
+
+global xap yap xdet ydet
+global ae da
+% ======== tokamak parameters ================================================
+
+load tcv_vesc1
+
+
+
+
+xcont=rzvin(:,1);
+ycont=rzvin(:,2);
+xmin=min(xcont);
+ymin=min(ycont);
+xmax=max(xcont);
+ymax=max(ycont);
+xedge=100;
+yedge=60;
+          
+
+
+
+% =========  detector parameters =============================================
+
+if i_detec==2
+
+	cw=1;                         % detector numbers cw=1:clockwise cw=0:ccw
+
+	if nargin<2
+		fans=[0 0 0 0 0 0 1 0 1 0];   % camera switch
+	end
+
+	vangle=[90 90 90 0 0 0 0 -90 -90 -90];
+                                % angle of detector surface normal 
+	xpos=[0 0 0 0 0 0  118.05 0 87.84 0]; 
+                                % x position of the diaphragmas in [cm]
+	ypos=[0 0 0 0 0 0 -46 0 -80.45 0];
+                                % y position of the diaphragmas in [cm]
+	ae=[0 0 0 0 0 0 -2.5 0 -0.1 0]/10; 
+				%excentricity array/diaphragm in [cm]
+	da=[0 0 0 0 0 0 10.1  0 -11.7 0]/10;   
+                                 % diaphragma-array distance in [cm]
+	da2=da;
+
+	d1=0.950;                     % detector width      in mm
+	d2=4.050;                     % detector length     in mm
+	b1=1.000;                     % aperture width      in mm
+	b2=4.000*ones(1,10);          % aperture length     in mm
+	b3x=0;                        % aperture thickness  in mm
+	b3y=0;
+
+elseif i_detec==1
+
+   cw=1;                         % detector numbers cw=1:clockwise cw=0:ccw
+
+   if nargin<2	
+   	fans=[1 1 1 1 1 1 1 1 0 1];   % camera switch
+   end
+
+   vangle=[90 90 90 0 0 0 0 -90 -90 -90];
+                                 % angle of detector surface normal 
+   xpos=[71.5 87.7 103.9 123.1  123.1 123.1 123.1 104.04 87.84 71.64]; 
+                                 % x position of the diaphragmas in [cm]
+   ypos=[80.35 80.35 80.35 48.1 1.95 -2.45 -48.6 -80.35 -80.35 -80.35];
+                                 % y position of the diaphragmas in [cm]
+   ae=[-8 0 8 5 9 -9 -5 8 0 -8]/10; %excentricity array/diaphragm [cm]
+
+
+  ae= ae + [-0.0915 0   0.1361    0.2123    0.0923   -0.0994   ...
+ 		 0.0872   -0.1520  0  0.9410 ]/10;
+  ae(1)=ae(1)+0.1/10;
+  ae(3)=8/10+0.14/10;
+  ae(4)=4.9/10;
+  ae(5)=9/10+0.2/10;
+  ae(6)=ae(6)-0.2/10;
+  ae(7)=-4.9/10;
+  ae(10)=-7.1/10;
+
+   da=   [12.4 9.9 12.4 9.9 13.4 13.4  9.9   -12.4  -9.9 -12.4]/10;
+                                 % diaphragma-array distance in [cm] (poloidal)
+   da2=[37 34.4 37 55.9 59.4 59.4 55.9 37 34.4 37]/10;
+                                 % dist to diaphragm in toroidal direction [cm];
+   deltada=[   -0.0311 0  -0.0458   -0.1179   -0.0615   -0.1105 ...
+			  -0.0510   -0.0515  0  -0.3223]/10;
+  deltada(4)=0;
+  deltada(6)=0;
+                                                                 
+
+   da=da+ deltada;
+
+   da2=da2+deltada;
+
+
+   d1=0.90;			 % detector width      in mm
+   d2=4.0;			 % detector length     in mm
+   b1=0.800;			 % aperture width      in mm (pol.)
+   b2=[8 8 8 15 15 15 15 8 8 8];			 
+				 % aperture length     in mm (tor.)
+   b3x=0.020;			 % aperture thickness  in mm (poloidal)
+   b3y=0;			 % aperture thickness  in mm (toroidal)
+end
+
+
+
+%======== calculation of the chords of view ===================================
+
+nact=sum(fans);
+iact=find(fans);
+ndet=20;
+ncam=10;
+
+
+
+% ---- apertures: ------------------
+
+xap=ones(ndet,1)*xpos(iact);
+xap=xap(:)';
+yap=ones(ndet,1)*ypos(iact);
+yap=yap(:)';
+
+% ---- detectors: ------------------
+
+vorz(find(vangle==90))=(-1)^(cw+1)*ones(size(find(vangle==90)));
+vorz(find(vangle==0))=(-1)^cw*ones(size(find(vangle==0)));
+vorz(find(vangle==-90))=(-1)^cw*ones(size(find(vangle==-90)));
+
+
+dete=(-9.025:0.950:9.025)'/10*vorz(iact)+ones(ndet,1)*ae(iact);
+
+dum_ae=dete(:)';
+
+dum_vangle=ones(ndet,1)*vangle(iact);
+dum_vangle=dum_vangle(:)';
+
+
+
+ivert=find(dum_vangle==90 | dum_vangle==-90);
+ihori=find(dum_vangle==0);
+
+dum_da=ones(ndet,1)*da(iact);
+dum_da=dum_da(:)';
+
+dxd=zeros(1,ndet*nact);
+dyd=zeros(1,ndet*nact);
+
+dxd(ivert)=dum_ae(ivert);
+dxd(ihori)=dum_da(ihori);
+
+dyd(ivert)=dum_da(ivert);
+dyd(ihori)=dum_ae(ihori);
+
+xdet=xap+dxd;
+ydet=yap+dyd;
+
+
+%plot_vessel(rzvin,rzvout)
+%hold on
+% plot(xap,yap,'.g',xdet,ydet,'.m')
+
+
+% ---- calculate the equations of lines of sight
+
+m=(ydet-yap)./(xdet-xap);
+b=ydet-m.*xdet;
+
+nl=length(xdet);
+xchord=zeros(2,nl);
+ychord=zeros(2,nl);
+
+
+xchord(1,:)=xdet;ychord(1,:)=ydet;
+
+iup=find(dum_vangle==90);
+isi=find(dum_vangle==0);
+ido=find(dum_vangle==-90);
+
+if ~isempty(iup)
+ychord(2,iup)=ymin*ones(size(iup));
+xchord(2,iup)=(ychord(2,iup)-b(iup))./m(iup);
+end
+if ~isempty(ido)
+ychord(2,ido)=ymax*ones(size(ido));
+xchord(2,ido)=(ychord(2,ido)-b(ido))./m(ido);
+end
+if ~isempty(isi)
+xchord(2,isi)=xmin*ones(size(isi));
+ychord(2,isi)=m(isi).*xchord(2,isi)+b(isi);
+end
+
+ileft=find(xchord(2,:)<xmin);
+
+if ~isempty(ileft)
+xchord(2,ileft)=xmin*ones(size(ileft));
+ychord(2,ileft)=m(ileft).*xchord(2,ileft)+b(ileft);
+end
+irig=find(xchord(2,:)>xmax);
+if ~isempty(irig)
+xchord(2,irig)=xmax*ones(size(irig));
+ychord(2,irig)=m(irig).*xchord(2,irig)+b(irig);
+end
+
+
+
+%======== prepare output ======================================================
+
+vangle=vangle(iact);
+
+%======== calculation of angular correction factors, if necessary =============
+
+if i_detec==2 & exist('angular_fact_2.mat')==2
+
+		disp('loading angular_fact_2')
+		load angular_fact_2
+
+
+elseif i_detec==1 & exist('angular_fact_1.mat')==2
+
+		disp('loading angular_fact_1')
+		load angular_fact_1
+
+else
+
+		aomega=zeros(ndet,ncam);
+		angfact=ones(ndet,ncam);
+
+
+		for l=1:sum(fans)
+
+%		Z0X=abs(da(iact(l))*10)
+%	 	Z0Y=abs(da2(iact(l))*10)
+%		X0=ae(iact(l))*10 % back to mm, sorry about that...
+%		X0=X0*vorz(iact(l))
+%		B2=b2(iact(l))
+%		AOMEGA=etendue_n(b1,B2,b3x,b3y,Z0X,Z0Y,X0,cw);
+
+		b1x=0.8;
+		b1y=6;
+		b1z=0.02;
+		b2x=10000000;
+		b2y=b2(iact(l));
+		b2z=0;
+		z01=abs(da(iact(l))*10);
+		z02=abs(da2(iact(l))*10);
+		X0=ae(iact(l))*10*vorz(iact(l));
+keyboard
+		AOMEGA=etendue_n2(b1x,b1y,b1z,b2x,b2y,b2z,z01,z02,X0,cw);
+
+			aomega(:,iact(l))=AOMEGA(:,1);
+
+		
+		end
+
+
+		indm=min(find(aomega==max(aomega(:))));
+		aomegan=aomega/aomega(indm);
+		nonz=find(aomega);
+		angfact(nonz)=ones(size(nonz))./aomegan(nonz);
+
+		angfact=round(1000*angfact)/1000;
+
+
+		unitstring='units aomega: mm^2 * sterad';
+
+		if i_detec==1
+			save angular_fact_1 angfact aomega unitstring
+		elseif i_detec==2
+			save angular_fact_2 angfact aomega unitstring
+		end
+
+end
+return
+
+th=atan(diff(ychord)./diff(xchord));
+thp=th;
+neg=find(thp<0);
+thp(neg)=180+thp(neg);
+thdet=ones(1,20*sum(fans));
+for k=1:sum(fans)
+     thdet((k-1)*20+1:k*20)=vangle(k)*ones(1,20);
+end
+angles=thdet-thp;
+mist=find(angles<0 & abs(angles)>90);
+angles(mist)=angles(mist)+180;
+th_inc=angles*pi/180;
+
+
+% ---- correct for the edges of tcv ( some chords may be too long )
+
+
+
+
+
+	down=find(xcont>xedge & ycont<-yedge);
+	up=find(xcont>xedge & ycont>yedge);
+	cd=polyfit(xcont(down),ycont(down),1);
+	cu=polyfit(xcont(up),ycont(up),1);
+
+
+	iu1=find(xchord(1,:)>xedge & ychord(1,:)>0 & dum_vangle==-90 );
+	if ~isempty(iu1)
+		xchord(1,iu1)=-(b(iu1)-cu(2))./(m(iu1)-cu(1)+eps);
+		ychord(1,iu1)=m(iu1).*xchord(1,iu1)+b(iu1);
+	end
+
+	iu2=find(xchord(2,:)>xedge & ychord(2,:)>0 & ychord(1,:) & ....
+			dum_vangle==-90);
+	if ~isempty(iu2)
+		xchord(2,iu2)=-(b(iu2)-cu(2))./(m(iu2)-cu(1)+eps);
+		ychord(2,iu2)=m(iu2).*xchord(2,iu2)+b(iu2);
+	end
+
+	id1=find(xchord(1,:)>xedge & ychord(1,:)<0 & dum_vangle==90);
+	if ~isempty(id1)
+		xchord(1,id1)=-(b(id1)-cd(2))./(m(id1)-cd(1)+eps);
+		ychord(1,id1)=m(id1).*xchord(1,id1)+b(id1);
+	end
+
+	id2=find(xchord(2,:)>xedge & ychord(2,:)<0 & dum_vangle==90);
+	if ~isempty(id2)
+		xchord(2,id2)=-(b(id2)-cd(2))./(m(id2)-cd(1)+eps);
+		ychord(2,id2)=m(id2).*xchord(2,id2)+b(id2);
+	end
+
+	ilow=find(ychord(1,:)<ymin);
+	ihig=find(ychord(1,:)>ymax);
+	ilef=find(xchord(1,:)<xmin);
+	irig=find(xchord(1,:)>xmax);
+	if ~isempty(ilow)
+		ychord(1,ilow)=ymin*ones(size(ilow));
+		xchord(1,ilow)=ymin./m(ilow)-b(ilow)./m(ilow);
+	end
+	if ~isempty(ihig)
+		ychord(1,ihig)=ymax*ones(size(ihig));
+		xchord(1,ihig)=ymax./m(ihig)-b(ihig)./m(ihig);
+	end
+	if ~isempty(ilef)
+		xchord(1,ilef)=xmin*ones(size(ilef));
+		ychord(1,ilef)=m(ilef)*xmin+b(ilef);
+	end
+	if ~isempty(irig)
+		xchord(1,irig)=xmax*ones(size(irig));
+		ychord(1,irig)=m(irig)*xmax+b(irig);
+	end
+
+
+	ilow=find(ychord(2,:)<ymin);
+	ihig=find(ychord(2,:)>ymax);
+	ilef=find(xchord(2,:)<xmin);
+	irig=find(xchord(2,:)>xmax);
+	if ~isempty(ilow)
+		ychord(2,ilow)=ymin*ones(size(ilow));
+		xchord(2,ilow)=ymin./m(ilow)-b(ilow)./m(ilow);
+	end
+	if ~isempty(ihig)
+		ychord(2,ihig)=ymax*ones(size(ihig));
+		xchord(2,ihig)=ymax./m(ihig)-b(ihig)./m(ihig);
+	end
+	if ~isempty(ilef)
+		xchord(2,ilef)=xmin*ones(size(ilef));
+		ychord(2,ilef)=m(ilef)*xmin+b(ilef);
+	end
+	if ~isempty(irig)
+		xchord(2,irig)=xmax*ones(size(irig));
+		ychord(2,irig)=m(irig)*xmax+b(irig);
+	end
+
+
+
+
+
-- 
GitLab