Commit e5e8a5fc authored by riccardo's avatar riccardo
Browse files

Lot of stuff

parent 3eb3632d
Loading
Loading
Loading
Loading
+29 −57
Original line number Diff line number Diff line
@@ -150,61 +150,26 @@ Omega_min = 2*pi/porb_max;
% singah = zeros(4,1);
nismin = zeros(4,1);
nismax = zeros(4,1);
% for s = 1:4
%     singah(s) = max(sin(gam + s*pi/2),[],'all');
%     singal(s) = min(sin(gam + s*pi/2),[],'all');
% end
% 
% %% Tutto sto macello e poi ora mi viene da pensare che valori intermedi di
% %% tasc possono portare a sin(gam) anche più alti o più bassi...
% %% Guarda come ora avrei fatto molto prima a mettere semplicemente +1 e -1
% % Adding special case for s = 1
% s = 1;
% if singal(s)>0
%     nismin(s) = -f_max*a_max*(Omega_max^s)*singah(s) + f_max;
%     nismax(s) = -f_min*a_min*(Omega_min^s)*singal(s) + f_min;
% elseif singah(s)>0
%     nismin(s) = -f_max*a_max*(Omega_max^s)*singah(s) + f_max;
%     nismax(s) = -f_max*a_max*(Omega_max^s)*singal(s) + f_max;
% else 
%     nismin(s) = -f_max*a_max*(Omega_max^s)*singah(s) + f_max;
%     nismax(s) = -f_min*a_min*(Omega_min^s)*singal(s) + f_min;
% end
% for s = 2:4
%     % This range is computed by finding the maximum span of Equation (15) after varying the search
%     % parameters over their respective ranges (given in Table 2). This
%     % is done with the exception of ν which is held fixed at its
%     % maximum value within sub-bands over the frequency search space.
%     %% RECHECK EVERY COMBINATION
%     if singal(s)>0
%         nismin(s) = -f_max*a_max*(Omega_max^s)*singah(s);
%         nismax(s) = -f_min*a_min*(Omega_min^s)*singal(s);
%     elseif singah(s)>0
%         nismin(s) = -f_max*a_max*(Omega_max^s)*singah(s);
%         nismax(s) = -f_max*a_max*(Omega_max^s)*singal(s);
%     else
%         nismin(s) = -f_max*a_max*(Omega_max^s)*singah(s);
%         nismax(s) = -f_min*a_min*(Omega_min^s)*singal(s);
%     end
% end

% Find the range in ν^s (i.e. the maximum span of Eq. (15) covered by 
% combining the search parameters over their respective ranges)
% Sin(gamma) going from -1 to 1
% Adding special case for s = 1
s = 1;
if a_max*Omega_max>1 
    nismin(s) = -f_max*a_max*(Omega_max) + f_max;
    if a_min*Omega_min>1
        nismax(s) = -f_min*a_min*(Omega_min) + f_min;
    else
        nismax(s) = -f_max*a_min*(Omega_min) + f_max;
    end
else
% if a_max*Omega_max>1 
%     nismin(s) = -f_max*a_max*(Omega_max) + f_max;
%     if a_min*Omega_min>1
%         nismax(s) = -f_min*a_min*(Omega_min) + f_min;
%     else
%         nismax(s) = -f_max*a_min*(Omega_min) + f_max;
%     end
% else
%     nismin(s) = -f_min*a_max*(Omega_max) + f_min;
%     nismax(s) = -f_max*a_min*(Omega_min) + f_max;
% end
% Is (1-Ome*a) more than fmin-fmax? to ask myself
nismin(s) = -f_min*a_max*(Omega_max) + f_min;
nismax(s) = -f_max*a_min*(Omega_min) + f_max;
end
% The other s's are easier
for s = 2:4
    nismin(s) = -f_max*a_max*(Omega_max^s);
    nismax(s) =  f_max*a_max*(Omega_max^s);
@@ -274,6 +239,12 @@ tic
tm=gpuArray(tm);
nibank = gpuArray(combinations(nis{:}).Variables);
Lambda = zeros(length(nibank),length(f_gr),'gpuArray');
F1=gpuArray((0:N-1)./(N*dt_psd)).';
f_ind = zeros(1,length(f_gr),'uint32');
for n = 1:length(f_gr)
    [ble,f_ind(n)] = min(abs(F1-f_gr(n)));
end
clear ble
% nibank = combinations(nis{:}).Variables;
% Lambda = zeros(length(nibank),length(f_gr));
toc 
@@ -304,7 +275,9 @@ for m=1:M
    ttemp = tm(m,:);
    ttempdif = (ttemp(:)-tmid(m)).';
    ttempdifl = (ttemp(:)-tmid(m)-0.5*dt).';
    xtemp = x(:,m).';
    xtemp = gpuArray(x(:,m).');
    sumx = sum(xtemp);
    

    % toc
    
@@ -338,23 +311,22 @@ for m=1:M
            taul(1:2) = tmid(m) + sum((nibank(i,1:s_s)./(nizero*factorial(1:s_s))).*((ttempdifl(1:2)).^((1:s_s).').'),2);
            dtau = ((taul(2)-taul(1))); 

            % inzo = (interp1(tau(:),(xtemp(:)./dtau),ttemp(:),'linear',0)).*dt;
            Y1 = fft((interp1(tau(:),(xtemp(:)./dtau),ttemp(:),'linear',0)).*dt).';
            F1=((0:length(Y1)-1)./(length(Y1)*dt_psd)).';
            [nnfreq,nnind] = min(abs(F1-nizero));
            Y1=Y1(nnind); clear F1
            Y1=Y1(f_ind(n)); 

            Lambda(i,n) = 2*(abs(max(Y1)).^2)/sum(xtemp(:)); %CREDO (oppure prendono la potenza massima?)
            Lambda(i,n) = 2*(abs(Y1)^2)/sumx; %CREDO (oppure prendono la potenza massima?)

            % toc
        end
        % disp('1ni')
    end
    toc
    disp(length(f_gr));
    disp(length(f_gr)*i);
    disp(m);
    lamfiname = sprintf('Lambda_seg_%d.mat', m);
    % save(lamfiname,"Lambda");
    save([pathfi,lamfiname],"Lambda");
    save([pathfi,lamfiname],"Lambda","-v7.3");
    Lambda = zeros(length(nibank),length(f_gr));
    disp(m);

@@ -409,7 +381,7 @@ toc

% disp(bestpar);
% save('C:\Users\Filippo\Desktop\XMM_Jxxx\risultelli.mat');
save([pathfi,'risultelli_',num2str(Tseg),'s.mat']);
save([pathfi,'risultelli_',num2str(Tseg),'s.mat'],"-v7.3");


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%