classdef (SharedTestFixtures={mexm_fixture}) locS_tests < meq_test
  % Tests of locSmex profile interpolation routines
  %
  % [+MEQ MatlabEQuilibrium Toolbox+] Swiss Plasma Center EPFL Lausanne 2022. All rights reserved.
  properties
    verbosity = 0;
  end
  
  properties(TestParameter)
    method = {'locSmex','locSmexm'};
    sq     = struct('positive',-1,'negative',1);
  end
  
  methods(Test,TestTags={'Unit'})
    
    function test_basic(testCase,method,sq)
      x = [0,0.5,1]; y = [0;0.5;1]; yy=[0.5;0.6;1];
      xx = feval(method,x,sq.*y,sq.*yy).';
      testCase.verifyEqual(yy,xx,'basic straight line interpolation failed')
    end
    
    function test_matrix(testCase,method,sq)
      X = [1;2]*[0,0.5,1]; y = [0;0.5;1]; yy=[0.5;0.6];

      XX = feval(method,X,sq.*y,sq.*yy); % direct call

      % same with for loop
      xx = zeros(size(XX));
      for io=1:size(X,1)
        x = X(io,:);
        xx(io,:)=feval(method,x,sq.*y,sq.*yy);
      end
      testCase.verifyEqual(xx,XX,'matrix and for loop don''t give same result');
    end
    
    function test_v_interp1(testCase,method,sq)
      %% A simple q profile
      aQ = 1 + linspace(0,1,21).^1.5; % unequally spaced grid
      qQ = 1 + 20*abs((aQ-1.4))'.^3; % some nonmonotonic profile
      
      % requested q values - some do not exist
      qS = [1,1.5,4,10];
      
      % amount of points to find
      
      % run locSmex() and variants
      aS = feval(method,aQ,sq./qQ,sq./qS);
      
      %% check vs interp1 for monotonic part
      imono = find(diff(qQ)>0)+1;
      aR1 = interp1(1./qQ(imono),aQ(imono),1./qS,'linear',NaN);
      
      testCase.verifyEqual(aS,aR1,'AbsTol',mean(diff(aQ))*1e-2)

      %% test empty
      aS = feval(method,aQ,sq./qQ,[]);
      testCase.assertEmpty(aS)
      
    end
  end
end
