Bfit_dke_ds

 Calculation of B and dBdtheta (option) on a psi grid for any poloidal
 angle value

 INPUT:

   - theta: poloidal angle [1,1]
   - psin: normalized poloidal flux grid [nr,1] 
   - equil: magnetic equilibrium structure
   - equil_fit: structure of the vectorial for of the magnetic equilibrium
   - Bsolve: reference magnetic field value (for use with thetastar_dke_ds.m only with fsolve matlab function) (T) [1,1]

 OUTPUT:

   - B_fit: magnetic field at theta value (T) [nr,1]
   - dBdtheta_fit: poloidal derivative of the magnetic field at theta value (T) [nr,1]
     (optional)

 by Y. PEYSSON (CEA/DRFC, yves.peysson@cea.fr) and D. SCHMITT (CEA/DRFC, damien.schmitt@cea.fr)

0001 function [varargout] = Bfit_dke_ds(theta,psin,equil,equil_fit,Bsolve)
0002 %
0003 % Calculation of B and dBdtheta (option) on a psi grid for any poloidal
0004 % angle value
0005 %
0006 % INPUT:
0007 %
0008 %   - theta: poloidal angle [1,1]
0009 %   - psin: normalized poloidal flux grid [nr,1]
0010 %   - equil: magnetic equilibrium structure
0011 %   - equil_fit: structure of the vectorial for of the magnetic equilibrium
0012 %   - Bsolve: reference magnetic field value (for use with thetastar_dke_ds.m only with fsolve matlab function) (T) [1,1]
0013 %
0014 % OUTPUT:
0015 %
0016 %   - B_fit: magnetic field at theta value (T) [nr,1]
0017 %   - dBdtheta_fit: poloidal derivative of the magnetic field at theta value (T) [nr,1]
0018 %     (optional)
0019 %
0020 % by Y. PEYSSON (CEA/DRFC, yves.peysson@cea.fr) and D. SCHMITT (CEA/DRFC, damien.schmitt@cea.fr)
0021 %
0022 if nargin < 5, Bsolve = 0;end
0023 %
0024 rho = psi2rho_jd(equil,psin);%psi to rho conversion
0025 %
0026 Bx_a0_fit = ppval_yp(equil_fit.Bx_fit.pp_a0,rho);
0027 Bx_an_fit = ppval_yp(equil_fit.Bx_fit.pp_an,rho);
0028 Bx_bn_fit = ppval_yp(equil_fit.Bx_fit.pp_bn,rho);
0029 %
0030 By_a0_fit = ppval_yp(equil_fit.By_fit.pp_a0,rho);
0031 By_an_fit = ppval_yp(equil_fit.By_fit.pp_an,rho);
0032 By_bn_fit = ppval_yp(equil_fit.By_fit.pp_bn,rho);
0033 %
0034 BPHI_a0_fit = ppval_yp(equil_fit.BPHI_fit.pp_a0,rho);
0035 BPHI_an_fit = ppval_yp(equil_fit.BPHI_fit.pp_an,rho);
0036 BPHI_bn_fit = ppval_yp(equil_fit.BPHI_fit.pp_bn,rho);
0037 %
0038 % Build interpolated quantities
0039 %
0040 [BR_fit,dBRdtheta_fit] = calcval_yp(equil_fit,theta,Bx_a0_fit(:),Bx_an_fit,Bx_bn_fit);
0041 [BZ_fit,dBZdtheta_fit] = calcval_yp(equil_fit,theta,By_a0_fit(:),By_an_fit,By_bn_fit);
0042 [BPHI_fit,dBPHIdtheta_fit] = calcval_yp(equil_fit,theta,BPHI_a0_fit(:),BPHI_an_fit,BPHI_bn_fit);
0043 %
0044 BP_fit = sqrt(BR_fit(:).^2 + BZ_fit(:).^2);
0045 dBPdtheta_fit = (BR_fit(:).*dBRdtheta_fit(:) + BZ_fit(:).*dBZdtheta_fit(:))./BP_fit(:);
0046 %
0047 B_fit = sqrt(BR_fit(:).^2 + BZ_fit(:).^2 + BPHI_fit(:).^2);
0048 dBdtheta_fit = (BR_fit(:).*dBRdtheta_fit(:) + BZ_fit(:).*dBZdtheta_fit(:) + BPHI_fit(:).*dBPHIdtheta_fit(:))./B_fit(:);
0049 %
0050 B_fitsolve = B_fit - Bsolve;%for finding theta value corresponding to B_fitsolve = 0 if Bsolve is defined.
0051 %
0052 if nargout >= 1, varargout{1} = B_fitsolve;end
0053 if nargout >= 2, varargout{2} = dBdtheta_fit;end