Commit ea84a252 authored by Riccardo La Placa's avatar Riccardo La Placa
Browse files

Pre-test

parent 652d5ddd
Loading
Loading
Loading
Loading
+29 −15
Original line number Diff line number Diff line
@@ -7,7 +7,8 @@ pathfi = './';
%%%%%%%%%%%%%%%%
tic

finame = 'J1023_B_2017_Bary.fits';
% finame = 'J1023_B_2017_Bary.fits';
finame = 'EPN_0744840201_bary.fits';
t_raw = fitsread([pathfi,finame],"binarytable");
t_raw = t_raw{1};
% info = fitsinfo(finame).BinaryTable.Keywords;
@@ -23,10 +24,17 @@ toc
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Lines added just to test what was done until now
f_tru = 592.42146827248556; %Hz
porb_tru = 17115.5216592; %s
a_tru = 0.343356; %lt-s
tasc_tru = 58107.009477; %MJD?
% Parameters for J1023_B_2017_Bary.fits
%f_tru = 592.42146827248556; %Hz
%porb_tru = 17115.5216592; %s
%a_tru = 0.343356; %lt-s
%tasc_tru = 58107.009477; %MJD?

% Parameters for EPN_0744840201_bary.fits
f_tru = 598.8921309; %Hz
porb_tru = 8844.08; %s
a_tru = 0.0649905; %lt-s
tasc_tru = 57231.437581; %MJD?

%f_gr = zeros(5,1);
porb_gr = zeros(5,1);
@@ -38,10 +46,10 @@ tic
f_gr=f_tru+(-2:2).';

for j = 1:5
    %f_gr(j) = f_tru*((j^3)/9);
    porb_gr(j) = porb_tru*((j^3)/9);
    a_gr(j) = a_tru*((j^3)/9);
    tasc_gr(j) = tasc_tru*((j^3)/9);
    %f_gr(j) = f_tru*((j^2)/9);
    porb_gr(j) = porb_tru*((j^2)/9);
    a_gr(j) = a_tru*((j^2)/9);
    tasc_gr(j) = tasc_tru*((j^2)/9);
end

f_min = min(f_gr);
@@ -212,7 +220,7 @@ toc
%for each segment (lavoro su tm)
for m=1:M

    tic
    % tic
 %   [Cm,edges]=(histcounts(x(m,:),round((tm(m,end)-tm(m,1))/dt_psd))); % R - convincitene
 %   edges=edges(end)-edges(1); %mi dà il tempo preciso di tutta la TdF, che sarà leggermente diversa da length(C)*dt per come è definito histcounts 
    % Y=fft(x(m,:)).'; 
@@ -225,7 +233,7 @@ for m=1:M
    % Y=Y(cond); 
    % X=ifft(Y); %inverse-fourier transf.

    toc
    % toc
    
    %FARE IL RICAMPIONAMENTO!!!!
    % Dev'essere fatto per ciascun template, ergo mettiamo un bel ciclo 
@@ -252,12 +260,14 @@ for m=1:M
        %         tau(j) = tau(j) + (nibank(i,s)/(nizero*factorial(s)))*(tm(m,j)-tmid(m))^s;
        %     end
        % end
        
        tic
      
        tau = sum((nibank(i,1:s_s)./(nizero*factorial(1:s_s))).*((tm(m,1:N)-tmid(m)).^((1:s_s).').'),2);
        toc

        tic
        X1 = interp1(tm(m,:),x(m,:),tau,'linear',0);
        
        toc

        %X1 è la timeseries ricampionata (controllare che sia un vettore colonna)
        %zero-padding (metto gli zeri alla fine) -------------------------------
@@ -268,11 +278,14 @@ for m=1:M
        %[Cm,edges]=(histcounts(X1,round((X1(end)-X1(1))/dt_psd)));
        %edges=edges(end)-edges(2); %mi dà il tempo preciso di tutta la TdF, che sarà leggermente diversa da length(C)*dt per come è definito histcounts
        % Y1=(2./sum(X1).*abs(fft(X1)).^2).'; %normalizzazione Leahy, giusto?????
        tic
        Y1 = fft(X1).';
        clear X1
        F1=((0:length(Y1)-1)./(tm(m,end)-tm(m,1))).';
        % F1=F1(1:round(length(F1)/2));
        % Y1=Y1(1:round(length(Y1)/2));
        toc
        tic
        cond = F1>=f_min & F1<=f_max;
        F1=F1(cond);
        Y1=Y1(cond);
@@ -280,7 +293,8 @@ for m=1:M
        %Calcolo della detection statistic
        %QUESTA SARA' DA INIZIALIZZARE ALL'ESTERNO tipo Lambda = zeros(length(template),1);
        Lambda(i,m)=sum(abs(Y1).^2)/sum(x(m,:)); %CREDO (oppure prendono la potenza massima?)
        
        toc
        disp('1ni')
    end
    toc