I wanted to rewrite my octave code to make it easier to read and make it a bit faster, just a simple vectorisation. I then run the two implementations and gets two very different results.
Did I fail my implementation?

clear; close all; clc;

xlist = 10:15;
ylist = (20:22);

% Making a LUT from x and y : LUT --> nConf | x | y
  count = 1;
for kx=1:length(xlist)
 for ky=1:length(ylist)
   LUT(count,:) = [count xlist(kx) ylist(ky)];
   count = count+1;
   end
end

% Making the data --> nConf | measurement
% data :  1 x1
%         2 x2
%        ...
%        15 x15
%         1 x16
%         2 x17
%        ...
Nrounds = 100;
data = zeros(Nrounds*length(LUT),2); % Let's have a multiple rounds/cycles measurements
data(:,2) = cos(1:length(data)); %Measurements
data(1:floor(Nrounds)*length(LUT),1) = repmat(LUT(:,1),[floor(Nrounds) 1]); % Fill the nConf column

% Parameter to compute a noise approximation
wlen = 3; % window length
opt = 0; %see movstd doc

% Debug parameter
debug_val = 2; % This allows to select either data(:,1) or data(:,2)

%% METHOD 1 : Reshaping the data in a 3D matrix with each column in the third axis is a conf with Nrounds layers
t1 = tic;
Mdata_rs = reshape(data(:,debug_val),[length(ylist) length(xlist) Nrounds]);
movSTDdata_rs = movstd(Mdata_rs,wlen,opt,3);
movSTDdata_rs = reshape(movSTDdata_rs,[length(ylist) length(xlist) Nrounds]); %movstd may swap dimensions
output1 = median(movSTDdata_rs,3); %Get a noise approximation ignoring the possible base line
toc(t1)

%% METHOD 2 : Doing the same thing but using find() for x and y values
t2 = tic;
output2 = zeros([length(ylist) length(xlist)])*nan;
for mx = 1:length(xlist)
  x_val = xlist(mx); % Fetch the research x value
  x_idx = find(LUT(:,2)==x_val); % Get at what rows it appears in the LUT
 for my = 1:length(ylist)
   y_val = ylist(my); % Fetch the research y value
   y_idx = find(LUT(x_idx,3)==y_val); % Get at what row both values appear in the LUT
   conf_val = LUT(x_idx(1)+y_idx-1); % Get the conf number
   data_idx = find(data(:,1)==conf_val); % Fetch all the measurement where this conf nuumber is used
   data_val = data(data_idx,debug_val); % Get the said measurements
   output2(my,mx) = median(movstd(data_val,wlen,opt,1)); % Do the same math
 end
end
toc(t2)

output1
output2

I exectuted the exact same file (not script, file) with Matlab, and got the expected result.
The result are the following :