0001 function data_pedestal = calc_pedestal_jd(equil,diary_mode)
0002
0003 if nargin < 2,
0004 diary_mode = false;
0005 end
0006 if nargin < 1,
0007 error('not enough arguments')
0008 elseif ischar(equil),
0009 load(equil);
0010 end
0011
0012
0013
0014
0015 prho_G = equil.ptx(:,1)/equil.ptx(end,1);
0016 prho_P = sqrt(equil.psi_apRp/equil.psi_apRp(end));
0017
0018 Temax = max(equil.pTe);
0019 nemax = max(equil.pne);
0020
0021 pTe = equil.pTe/Temax;
0022 pne = equil.pne/nemax;
0023
0024 np = length(prho_G);
0025
0026 i_LFS = 1;
0027 i_HFS = (length(equil.theta) + 1)/2;
0028
0029 pB_LFS = sqrt(equil.ptBx(:,i_LFS).^2 + equil.ptBy(:,i_LFS).^2 + equil.ptBPHI(:,i_LFS).^2);
0030 pB_HFS = sqrt(equil.ptBx(:,i_HFS).^2 + equil.ptBy(:,i_HFS).^2 + equil.ptBPHI(:,i_HFS).^2);
0031
0032 pdnedrho = NaN(1,np);
0033 pdTedrho = NaN(1,np);
0034
0035 drhop = prho_P(3:end) - prho_P(2:end-1);
0036 drhom = prho_P(2:end-1) - prho_P(1:end-2);
0037
0038 dnep = pne(3:end) - pne(2:end-1);
0039 dnem = pne(2:end-1) - pne(1:end-2);
0040
0041 dTep = pTe(3:end) - pTe(2:end-1);
0042 dTem = pTe(2:end-1) - pTe(1:end-2);
0043
0044 pdnedrho(2:end - 1) = (drhop.^2.*dnem + drhom.^2.*dnep)./drhom./drhop./(drhop + drhom);
0045 pdnedrho(1) = 0;
0046 pdnedrho(end) = dnep(end)/drhop(end);
0047
0048 pdTedrho(2:end - 1) = (drhop.^2.*dTem + drhom.^2.*dTep)./drhom./drhop./(drhop + drhom);
0049 pdTedrho(1) = 0;
0050 pdTedrho(end) = dTep(end)/drhop(end);
0051
0052 figure(1),clf
0053 set(1,'name','profiles')
0054
0055 xlab = 'r/a';
0056 ylab = '';
0057 xlim1 = [0,1];
0058 ylim1 = [0,ceil(max([nemax/1e19,Temax]))+1];
0059 tit = '';
0060 leg = {'ne19','Te (keV)'};
0061
0062 graph1D_jd(prho_P,pne*nemax/1e19,0,0,'','','',NaN,xlim1,ylim1,'-','+','b',2,20+14i);
0063 graph1D_jd(prho_P,pTe*Temax,0,0,xlab,ylab,tit,leg,xlim1,ylim1,'-','s','r',2,20+14i,gca,0.9,0.7,0.7);
0064
0065 figure(2),clf
0066 set(2,'name','profiles and derivatives')
0067
0068 xlim2 = [0.75,1];
0069 ylim2 = [-0.1,1];
0070 tit = 'Select Pedestal location';
0071 leg = {'ne','Te','-dne/drho (x 0.1)','-dTe/drho (x 0.1)'};
0072
0073 graph1D_jd(prho_P,pne,0,0,'','','',NaN,xlim2,ylim2,'-','+','b',2,20+14i);
0074 graph1D_jd(prho_P,pTe,0,0,'','','',NaN,xlim2,ylim2,'-','+','r',2,20+14i);
0075 graph1D_jd(prho_P,-pdnedrho/10,0,0,'','','',NaN,xlim2,ylim2,'-','+',[0,0.5,0],2,20+14i);
0076 graph1D_jd(prho_P,-pdTedrho/10,0,0,xlab,ylab,tit,leg,xlim2,ylim2,'-','+','m',2,20+14i,gca,0.9,0.7,0.5);
0077
0078 rho = ginput(2);
0079
0080 irho = find(abs(prho_P - rho(1)) == min(abs(prho_P - rho(1))),1,'first');
0081
0082 figure(1),
0083 graph1D_jd(prho_P(irho)*[1,1],ylim1,0,0,'','','',NaN,xlim1,ylim1,'-','none','k',2,20+14i);
0084
0085 figure(2),
0086 graph1D_jd(prho_P(irho)*[1,1],ylim2,0,0,'','','',NaN,xlim2,ylim2,'-','none','k',2,20+14i);
0087
0088 if diary_mode
0089 print_jd(2,['Fig_',equil.id,'_pedestal_1'],NaN,1)
0090 print_jd(2,['Fig_',equil.id,'_pedestal_2'],NaN,2)
0091
0092 filename = ['PROFILES_EQUIL_',equil.id];
0093
0094 if exist([filename,'.txt'],'file'),
0095 delete([filename,'.txt'])
0096 end
0097 diary([filename,'.txt']);
0098 end
0099
0100 data_pedestal.position.rho_G = prho_G(irho);
0101 data_pedestal.position.rho_P = prho_P(irho);
0102
0103 data_pedestal.bhoriz.edge_HFS = pB_HFS(end);
0104 data_pedestal.bhoriz.ped_HFS = pB_HFS(irho);
0105 data_pedestal.bhoriz.core = pB_LFS(1);
0106 data_pedestal.bhoriz.ped_LFS = pB_LFS(irho);
0107 data_pedestal.bhoriz.edge_LFS = pB_LFS(end);
0108
0109 data_pedestal.ne20.max = nemax/1e20;
0110 data_pedestal.ne20.core = pne(1)*nemax/1e20;
0111 data_pedestal.ne20.ped = pne(irho)*nemax/1e20;
0112 data_pedestal.ne20.edge = pne(end)*nemax/1e20;
0113
0114 data_pedestal.tekeV.max = Temax;
0115 data_pedestal.tekeV.core = pTe(1)*Temax;
0116 data_pedestal.tekeV.ped = pTe(irho)*Temax;
0117 data_pedestal.tekeV.edge = pTe(end)*Temax;
0118
0119 disp('-------------------------')
0120 disp(' PROFILES PARAMETERS ')
0121 disp('-------------------------')
0122 disp(' ');
0123 disp('Pedestal : ');
0124 disp([' - rho G : ',num2str(data_pedestal.position.rho_G)]);
0125 disp([' - rho P : ',num2str(data_pedestal.position.rho_P)]);
0126 disp(' ');
0127 disp('Density (ne20) : ')
0128 disp([' - Maximum : ',num2str(data_pedestal.ne20.max)]);
0129 disp([' - Central : ',num2str(data_pedestal.ne20.core)]);
0130 disp([' - Pedestal : ',num2str(data_pedestal.ne20.ped)]);
0131 disp([' - Edge : ',num2str(data_pedestal.ne20.edge)]);
0132 disp(' ');
0133 disp('Temperature (keV) : ')
0134 disp([' - Maximum : ',num2str(data_pedestal.tekeV.max)]);
0135 disp([' - Central : ',num2str(data_pedestal.tekeV.core)]);
0136 disp([' - Pedestal : ',num2str(data_pedestal.tekeV.ped)]);
0137 disp([' - Edge : ',num2str(data_pedestal.tekeV.edge)]);
0138 disp(' ');
0139 disp('Magnetic field (T) : ')
0140 disp([' - Edge (HFS) : ',num2str(data_pedestal.bhoriz.edge_HFS)]);
0141 disp([' - Pedestal (HFS) : ',num2str(data_pedestal.bhoriz.ped_HFS)]);
0142 disp([' - Central : ',num2str(data_pedestal.bhoriz.core)]);
0143 disp([' - Pedestal (LFS) : ',num2str(data_pedestal.bhoriz.ped_LFS)]);
0144 disp([' - Edge (LFS) : ',num2str(data_pedestal.bhoriz.edge_LFS)]);
0145
0146 disp(' ');
0147 disp('-------------------------')
0148
0149 if diary_mode
0150 diary off
0151
0152 save([filename,'.mat'],'data_pedestal');
0153 end
0154
0155
0156
0157
0158
0159
0160
0161
0162