0001 function equil_RZ = pt2rz_jd(equil,x,y,nprz,method)
0002
0003 [np,nt] = size(equil.ptx);
0004
0005 if nargin < 5,
0006 method = 'linear';
0007 end
0008
0009 if nargin < 4 || isempty(nprz),
0010 nprz = np;
0011 end
0012
0013 equil_RZ.psi_apRp = equil.psi_apRp;
0014 equil_RZ.ap = equil.ptx(end,1);
0015 equil_RZ.Rp = equil.Rp;
0016 equil_RZ.Zp = equil.Zp;
0017 if isfield(equil,'zZi')
0018 equil_RZ.zZi = equil.zZi;
0019 end
0020 if isfield(equil,'zmi')
0021 equil_RZ.zmi = equil.zmi;
0022 end
0023 if isfield(equil,'id')
0024 equil_RZ.id = equil.id;
0025 else
0026 equil.id = '';
0027 equil_RZ.id = '';
0028 end
0029
0030 ptx = equil.ptx;
0031 pty = equil.pty;
0032
0033 if nargin < 3,
0034 y = [min(pty(end,:)),max(pty(end,:))];
0035 end
0036
0037 if nargin < 2,
0038 x = [min(ptx(end,:)),max(ptx(end,:))];
0039 end
0040
0041 if length(x) == 1,
0042 extfacx = x;
0043 x = NaN;
0044 end
0045
0046 if length(y) == 1,
0047 extfacy = y;
0048 y = NaN;
0049 end
0050
0051 psia_apRp = equil.psi_apRp(end);
0052 psin = equil.psi_apRp/psia_apRp;
0053
0054 dpsin = psin(end) - psin(end-1);
0055
0056 if isnan(x),
0057 xmin = min(ptx(end,:));
0058 xmax = max(ptx(end,:));
0059
0060 itxmin = find(ptx(end,:) == xmin,1,'first');
0061 itxmax = find(ptx(end,:) == xmax,1,'first');
0062
0063 dxmin = (ptx(end,itxmin) - ptx(end-1,itxmin))*extfacx/dpsin;
0064 dxmax = (ptx(end,itxmax) - ptx(end-1,itxmax))*extfacx/dpsin;
0065
0066 x = [xmin + dxmin,xmax + dxmax];
0067 end
0068
0069 if isnan(y),
0070 ymin = min(pty(end,:));
0071 ymax = max(pty(end,:));
0072
0073 itymin = find(pty(end,:) == ymin,1,'first');
0074 itymax = find(pty(end,:) == ymax,1,'first');
0075
0076 dymin = (pty(end,itymin) - pty(end-1,itymin))*extfacy/dpsin;
0077 dymax = (pty(end,itymax) - pty(end-1,itymax))*extfacy/dpsin;
0078
0079 y = [ymin + dymin,ymax + dymax];
0080 end
0081
0082 if min(x) < -equil.Rp,
0083 error('min(x) < -equil.Rp');
0084 end
0085
0086 dsdpsin_min = min(sqrt((ptx(end,:) - ptx(end-1,:)).^2 + (pty(end,:) - pty(end-1,:)).^2))/dpsin;
0087 psin_max = 1 + sqrt((max(x) - min(x)).^2 + (max(y) - min(y)).^2)/dsdpsin_min;
0088 psin_ext = [psin,1 + dpsin:dpsin:psin_max];
0089
0090 if length(x) > 2,
0091 equil_RZ.x = x;
0092 else
0093 equil_RZ.x = linspace(x(1),x(2),nprz);
0094 end
0095
0096 if length(y) > 2,
0097 equil_RZ.y = y;
0098 else
0099 equil_RZ.y = linspace(y(1),y(2),nprz);
0100 end
0101
0102 ptx_ext = interp1(psin,ptx,psin_ext,'linear','extrap');
0103 pty_ext = interp1(psin,pty,psin_ext,'linear','extrap');
0104 ptBx_ext = interp1(psin,equil.ptBx,psin_ext,'linear','extrap');
0105 ptBy_ext = interp1(psin,equil.ptBy,psin_ext,'linear','extrap');
0106 ptBPHI_ext = interp1(psin,equil.ptBPHI,psin_ext,'linear','extrap');
0107
0108 if isfield(equil,'ptr')
0109 ptr_ext = interp1(psin,equil.ptr,psin_ext,'linear','extrap');
0110 end
0111
0112 if isfield(equil,'ptB')
0113 ptB_ext = interp1(psin,equil.ptB,psin_ext,'linear','extrap');
0114 end
0115
0116 if isfield(equil,'ptdBdtheta')
0117 ptdBdtheta_ext = interp1(psin,equil.ptdBdtheta,psin_ext,'linear','extrap');
0118 end
0119
0120 if isfield(equil,'ptBP')
0121 ptBP_ext = interp1(psin,equil.ptBP,psin_ext,'linear','extrap');
0122 end
0123
0124 if isfield(equil,'ptdBPdtheta')
0125 ptdBPdtheta_ext = interp1(psin,equil.ptdBPdtheta,psin_ext,'linear','extrap');
0126 end
0127
0128
0129 if isfield(equil,'ptcalpha')
0130 ptcalpha_ext = interp1(psin,equil.ptcalpha,psin_ext,'linear','extrap');
0131 end
0132
0133 if isfield(equil,'ptsalpha')
0134 ptsalpha_ext = interp1(psin,equil.ptsalpha,psin_ext,'linear','extrap');
0135 end
0136
0137 if isfield(equil,'ptgradrho')
0138 ptgradrho_ext = interp1(psin,equil.ptgradrho,psin_ext,'linear','extrap');
0139 end
0140
0141 if isfield(equil,'ptUpsilon')
0142 ptUpsilon_ext = interp1(psin,equil.ptUpsilon,psin_ext,'linear','extrap');
0143 end
0144
0145 psi_apRp_ext = psia_apRp*repmat(psin_ext.',[1,nt]);
0146
0147 warning('off')
0148
0149 equil_RZ.Bx = griddata(ptx_ext,pty_ext,ptBx_ext,equil_RZ.x,equil_RZ.y.',method).';
0150 equil_RZ.By = griddata(ptx_ext,pty_ext,ptBy_ext,equil_RZ.x,equil_RZ.y.',method).';
0151 equil_RZ.BPHI = griddata(ptx_ext,pty_ext,ptBPHI_ext,equil_RZ.x,equil_RZ.y.',method).';
0152 equil_RZ.xypsi_apRp = griddata(ptx_ext,pty_ext,psi_apRp_ext,equil_RZ.x,equil_RZ.y.',method).';
0153
0154 if isfield(equil,'ptr')
0155 equil_RZ.r = griddata(ptx_ext,pty_ext,ptr_ext,equil_RZ.x,equil_RZ.y.',method).';
0156 end
0157
0158 if isfield(equil,'ptB')
0159 equil_RZ.B = griddata(ptx_ext,pty_ext,ptB_ext,equil_RZ.x,equil_RZ.y.',method).';
0160 end
0161
0162 if isfield(equil,'ptdBdtheta')
0163 equil_RZ.dBdtheta = griddata(ptx_ext,pty_ext,ptdBdtheta_ext,equil_RZ.x,equil_RZ.y.',method).';
0164 end
0165
0166 if isfield(equil,'ptBP')
0167 equil_RZ.BP = griddata(ptx_ext,pty_ext,ptBP_ext,equil_RZ.x,equil_RZ.y.',method).';
0168 end
0169
0170 if isfield(equil,'ptdBPdtheta')
0171 equil_RZ.dBPdtheta = griddata(ptx_ext,pty_ext,ptdBPdtheta_ext,equil_RZ.x,equil_RZ.y.',method).';
0172 end
0173
0174 if isfield(equil,'ptcalpha')
0175 equil_RZ.calpha = griddata(ptx_ext,pty_ext,ptcalpha_ext,equil_RZ.x,equil_RZ.y.',method).';
0176 end
0177
0178 if isfield(equil,'ptsalpha')
0179 equil_RZ.salpha = griddata(ptx_ext,pty_ext,ptsalpha_ext,equil_RZ.x,equil_RZ.y.',method).';
0180 end
0181
0182 if isfield(equil,'ptgradrho')
0183 equil_RZ.gradrho = griddata(ptx_ext,pty_ext,ptgradrho_ext,equil_RZ.x,equil_RZ.y.',method).';
0184 end
0185
0186 if isfield(equil,'ptUpsilon')
0187 equil_RZ.Upsilon = griddata(ptx_ext,pty_ext,ptUpsilon_ext,equil_RZ.x,equil_RZ.y.',method).';
0188 end
0189
0190 warning('on')
0191
0192 if isfield(equil,'pTe')
0193 equil_RZ.pTe = interp1(psin,equil.pTe,psin_ext,'linear','extrap');
0194 end
0195 if isfield(equil,'pne')
0196 equil_RZ.pne = interp1(psin,equil.pne,psin_ext,'linear','extrap');
0197 end
0198 if isfield(equil,'pzTi')
0199 equil_RZ.pzTi = interp1(psin,equil.pzTi.',psin_ext,'linear','extrap').';
0200 end
0201 if isfield(equil,'pzni')
0202 equil_RZ.pzni = interp1(psin,equil.pzni.',psin_ext,'linear','extrap').';
0203 end
0204
0205 if isfield(equil,'pTe')
0206 equil_RZ.pTe(equil_RZ.pTe < 0) = 0;
0207 end
0208 if isfield(equil,'pne')
0209 equil_RZ.pne(equil_RZ.pne < 0) = 0;
0210 end
0211 if isfield(equil,'pzTi')
0212 equil_RZ.pzTi(equil_RZ.pzTi < 0) = 0;
0213 end
0214 if isfield(equil,'pzni')
0215 equil_RZ.pzni(equil_RZ.pzni < 0) = 0;
0216 end
0217
0218 if isfield(equil,'pTe') & isfield(equil,'pne') & isfield(equil,'pzTi') & isfield(equil,'pzni'),
0219 mask = equil_RZ.pTe == 0 & equil_RZ.pne == 0 & sum(equil_RZ.pzTi,1) == 0 & sum(equil_RZ.pzni,1) == 0;
0220 psin_ext(mask) = [];
0221 end
0222
0223 if isfield(equil,'pTe')
0224 equil_RZ.pTe(mask) = [];
0225 end
0226 if isfield(equil,'pne')
0227 equil_RZ.pne(mask) = [];
0228 end
0229 if isfield(equil,'pzTi')
0230 equil_RZ.pzTi(:,mask) = [];
0231 end
0232 if isfield(equil,'pzni')
0233 equil_RZ.pzni(:,mask) = [];
0234 end
0235 equil_RZ.psi_apRp = psin_ext*psia_apRp;
0236
0237
0238
0239
0240