%% Funzione che calcola gli spostamenti in modalità biassiale per i 
% Tilt Link HR (ampolle)

function [X_IPLHR,Y_IPLHR,Z_IPLHR,Xlocal_IPLHR,Ylocal_IPLHR,Zlocal_IPLHR,...
    HShift_IPLHR,HShift_local_IPLHR,AlfaX_IPLHR,AlfaY_IPLHR,Azimuth_IPLHR,...
    Speed_local_IPLHR,Speed_IPLHR,Acceleration_local_IPLHR,Acceleration_IPLHR,...
    TempDef_IPLHR,ARRAYdateIPLHR,ErrInPlaceLinkHR] = biax_IPLHR(IDcentralina,...
    DTcatena,rIPLHR,rIPL,ANGdef_IPLHR,ACCdef_IPL,MAGdef_IPL,TempDef_IPLHR,...
    SpeIPLHR,PsIPLHR,yesIPL,yesIPLHR3D,NodoInPlaceLinkHR,NodoInPlaceLink,...
    NodoInPlaceLinkHR3D,DatiElabInPlaceLinkHR,segnoNS_HR,segnoEO_HR,Ndevst_HR,...
    Wdevst_HR,ARRAYdateIPLHR,NuovoZeroIPLHR,NdatiMedia,Ndatidespike,allineato,...
    Tmax,Tmin,datainiIPLHR,ErrInPlaceLinkHR,Corr_Azimuth,margine,FileName)

%% Inizializzazione
fileID = fopen(FileName,'a');
fmt = '%s \r';
text = 'biax_IPLHR function started';
fprintf(fileID,fmt,text);

if NuovoZeroIPLHR == 1
   if NdatiMedia > Ndatidespike
        Ndati = NdatiMedia;
    else
        Ndati = Ndatidespike;
    end
    ini = round(Ndati/2)+1;
    if rem(Ndati,2) == 0
        ini = ini+1;
    end
    clear NDati
    ini = ini + margine;
    if ini < 6
        ini = 6;
    end
    if Ndevst_HR ~= 0 % Allora prendo tutti i dati e solo in seguito considero i nuovi, a valle della funzione filtro
        ini = 1;
    end
    if strcmp(NodoInPlaceLinkHR(1,4),'ADC') || strcmp(NodoInPlaceLinkHR{1,4} ,'null') || isempty(NodoInPlaceLinkHR{1,4}) == 1
        [rA,~] = size(ANGdef_IPLHR); % Numero di dati delle ampolle
        [rM,~] = size(ACCdef_IPL); % Numero di dati dei MEMS
        if rA~=rM
            if rM > rA
                ini_MEMS = rM-rA+ini;
                ini_Ampolle = ini;
            elseif rM < rA
                ini_Ampolle = rA-rM+ini; % Si presume che la data finale sia identica per entrambi
                ini_MEMS = ini;
            end
        else
            ini_MEMS = ini;
            ini_Ampolle = ini;
        end
        MAGdef_IPL = MAGdef_IPL(ini_MEMS:end,:);
    elseif strcmp(NodoInPlaceLinkHR{1,4} ,'digit')
        ini_Ampolle = ini;
    end
    ANGdef_IPLHR = ANGdef_IPLHR(ini_Ampolle:end,:);
    TempDef_IPLHR = TempDef_IPLHR(ini_Ampolle:end,:);
    DatiElabInPlaceLinkHR = DatiElabInPlaceLinkHR(ini_Ampolle:end,:);
    ARRAYdateIPLHR = ARRAYdateIPLHR(ini_Ampolle:end,1); 
    ErrInPlaceLinkHR = ErrInPlaceLinkHR(ini_Ampolle:end,:);
end

%% Elaborazione sfruttando il magnetometro del MEMS
if allineato == 0 && yesIPLHR3D == 1   
    % Definisco i dati MEMS
    [r,~] = size(MAGdef_IPL);
    Nnodi = rIPL;
    mx = zeros(r,Nnodi);
    my = zeros(r,Nnodi);
    for i=1:Nnodi
        mx(:,i)= -1*MAGdef_IPL(:,(i-1)*3+1); % mx
        my(:,i) = -1*MAGdef_IPL(:,(i-1)*3+2:(i-1)*3+2); % my
    end

    mx = mx'; % riga nodi, colonna date
    my = my'; % riga nodi, colonna date
    
    NodoInPlaceLink = cell2mat(NodoInPlaceLink(:,2));
    NodoDoppio = cell2mat(NodoInPlaceLinkHR3D(:,2));
    [S] = size(NodoDoppio);
    Nodes = zeros(S(1),1);
    for s = 1:S(1)
        Nodes(s,1) = find(NodoInPlaceLink == NodoDoppio(s,1)); 
    end
end

% Definisco i dati della cella elettrolitica
[Ndate,~] = size(ARRAYdateIPLHR);
Nnodi = rIPLHR;
AngoloX = zeros(Ndate,Nnodi); % angoli di inclinazione asse X
AngoloY = zeros(Ndate,Nnodi); % angoli di inclinazione asse Y
ii = 1;
for i=1:Nnodi
    AngoloX(:,i) = ANGdef_IPLHR(:,ii); % ax
    AngoloY(:,i) = ANGdef_IPLHR(:,ii+1);% ay
    ii = ii+2;
end

AngoloX = AngoloX';
AngoloY = AngoloY';

%% Costruzione delle matrici spostamento e rotazione
sX = zeros(rIPLHR,Ndate); % in riga i nodi, in colonna le date
sY = zeros(rIPLHR,Ndate);
sZ = zeros(rIPLHR,Ndate);

% Inizio del ciclo di elaborazione
text = 'Elaboration of In Place Link HR started';
fprintf(fileID,fmt,text);

for Imis = 1:Ndate % date
    for jjj=1:Nnodi % nodi
        %% Uso magnetometro MEMS
        if allineato == 0 && yesIPLHR3D == 1 % uso il magnetometro del MEMS
            Mx = mx(Nodes(jjj,1),:);
            My = my(Nodes(jjj,1),:);
            angle = arotHR(Mx,My,Imis);
            % angolo di inclinazione rispetto all'orizzontale registrato da asse x
            xAngTLHR = AngoloX(jjj,Imis);
            % angolo di inclinazione rispetto all'orizzontale registrato da asse y
            yAngTLHR = AngoloY(jjj,Imis);
            % converto gli angoli in radianti (matlab ragiona così)
            xAngTLHR = 0.01745329251994329576923690768489*xAngTLHR;
            yAngTLHR = 0.01745329251994329576923690768489*yAngTLHR;
            % Segmento di pertinenza
            if yesIPLHR3D == 1
                SpIPLHR = SpeIPLHR(jjj+1);  % correzione che salta l'ancora
            else
                % Correzione che salta l'ancora e il primo nodo Tilt Link (se
                % presente)
                SpIPLHR = SpeIPLHR(jjj+2);
            end
            % Finalizzo
            [Na,Ea,Za] = ASSEa_HR(xAngTLHR,angle,SpIPLHR);
            [Nb,Eb,Zb] = ASSEb_HR(yAngTLHR,angle,SpIPLHR);
            sX(jjj,Imis) = Na+Nb;
            sY(jjj,Imis) = Ea+Eb;
            % Calcolo l'abbassamento (positivo) considerando il MAGGIORE dei due calcolati sui due assi
            sZ(jjj,Imis) = max(Za,Zb);
        else
            %% NON uso magnetometro MEMS o IPTM
            
            % angolo di inclinazione rispetto all'orizzontale registrato da asse x
            xAngTLHR = AngoloX(jjj,Imis);
            % angolo di inclinazione rispetto all'orizzontale registrato da asse y
            yAngTLHR = AngoloY(jjj,Imis);
            % converto gli angoli in radianti (matlab ragiona così)
            xAngTLHR = 0.01745329251994329576923690768489*xAngTLHR;
            yAngTLHR = 0.01745329251994329576923690768489*yAngTLHR;
            % Segmento di pertinenza
            if  yesIPL == 1
                if yesIPLHR3D == 1
                    SpIPLHR = SpeIPLHR(jjj+1);  % correzione che salta l'ancora
                else
                % Correzione che salta l'ancora e il primo nodo Tilt Link (se
                % presente)
                    SpIPLHR = SpeIPLHR(jjj+2); 
                end
            else
                SpIPLHR = SpeIPLHR(jjj+1);  % correzione che salta l'ancora
            end    
            % calcolo spostamenti lungo x e y
            sX(jjj,Imis) = SpIPLHR*sin(xAngTLHR); 
            sY(jjj,Imis) = SpIPLHR*sin(yAngTLHR);            
            % Abbassamento registrato dai due assi
            sXZ = SpIPLHR - SpIPLHR*cos(xAngTLHR); 
            sYZ = SpIPLHR - SpIPLHR*cos(yAngTLHR);
            % Calcolo l'abbassamento considerando il maggiore dei due calcolati sui due assi
            sZ(jjj,Imis) = max(sXZ,sYZ);
        end
    end
end

% dsX, dsY e dsZ raccolgono i dati del singolo nodo nella singola data
dsX = diff(sX,1,2);
dsY = diff(sY,1,2);
dsZ = diff(sZ,1,2);

%% Cambio di segno le direzioni! 
if segnoNS_HR == 0
    dsX= -1*dsX;
end
if segnoEO_HR == 0
    dsY = -1*dsY;
end

%% Filtro su temperatura
clear i
clear j
cont2 = 1; % contatore
TempDef_IPLHR = TempDef_IPLHR';
textT = 'There are not correction of In Place Link HR based on temperature filter';

FileTemperature = ['' IDcentralina '-' DTcatena '-IPLHR-Therm.csv'];
if isfile(FileTemperature) == 1
    DatiRaw = csvread(FileTemperature);
    [rDR,cDR] = size(DatiRaw);
    DatiRaw(:,1) = DatiRaw(:,1) + 730000;
else
    rDR = 1;
    cDR = 1;
end
for j = 1:Ndate % Data
    for i = 1:Nnodi % Nodo
        % NON considero i dati al di sopra di Tmax o al di sotto di Tmin
        if TempDef_IPLHR(i,j) > Tmax || TempDef_IPLHR(i,j) < Tmin
            cont2 = cont2+1;
            if j == 1
                if isfile(FileTemperature) == 1
                    RawDate = find(DatiRaw(:,1)<=datenum(datainiIPLHR));
                    if isempty(RawDate) == 1
                        cc = 2;
                        while cc <= Ndate
                            if TempDef_IPLHR(i,cc) > Tmax || TempDef_IPLHR(i,cc) < Tmin
                                cc = cc+1;
                            else
                                break
                            end
                        end
                        TempDef_IPLHR(i,j) = TempDef_IPLHR(i,cc);
                    else
                        if isnan(DatiRaw(RawDate(end),i+1)) == 0
                            TempDef_IPLHR(i,j) = DatiRaw(RawDate(end),i+1);
                            ErrInPlaceLinkHR(j,i) = 0.5;
                            wardat = 'Temperature data of In Place Link HR nodes corrected using Raw Data of reference Csv file.';
                            fprintf(fileID,fmt,wardat);
                        else
                            cc = 2;
                            while cc <= c
                                if TempDef_IPLHR(i,cc) > Tmax || TempDef_IPLHR(i,cc) < Tmin
                                    cc = cc+1;
                                else
                                    break
                                end
                            end
                            TempDef_IPLHR(i,j) = TempDef_IPLHR(i,cc);
                        end
                    end
                else
                    cc = 2;
                    while cc <= c
                        if TempDef_IPLHR(i,cc) > Tmax || TempDef_IPLHR(i,cc) < Tmin
                            cc = cc+1;
                        else
                            break
                        end
                    end
                    TempDef_IPLHR(i,j) = TempDef_IPLHR(i,cc);
                end
            else       
                dsX(i,j-1) = 0;
                dsY(i,j-1) = 0;
                dsZ(i,j-1) = 0;
                TempDef_IPLHR(i,j) = TempDef_IPLHR(i,j-1);
                ErrInPlaceLinkHR(j,i) = 0.5;
            end
        end
        textT = ['' num2str(cont2) ' correction executed for In Place Link HR - Temperature filter!'];
    end
end

if  rDR~=1 && cDR~=1 && isempty(DatiRaw) == 0
    RawDate1 = find(DatiRaw(:,1)<=ARRAYdateIPLHR(1));
    if isempty(RawDate1) == 1
        RawDate2 = 1;
    elseif RawDate1(end) == rDR
        RawDate2 = find(ARRAYdateIPLHR(:,1)>DatiRaw(end,1));
    else
        RawDate2 = find(ARRAYdateIPLHR(:,1)>DatiRaw(RawDate1(end)+1,1));
    end
else
    RawDate1 = [];
    RawDate2 = 1;
end
if isempty(RawDate1) == 0 && isempty(RawDate2) == 0
    Dati = [DatiRaw(1:RawDate1(end),:); ARRAYdateIPLHR(RawDate2(1):end) TempDef_IPLHR(:,RawDate2(1):end)'];
elseif isempty(RawDate1) == 1 && isempty(RawDate2) == 0
    Dati = [ARRAYdateIPLHR TempDef_IPLHR'];
else
    Dati = DatiRaw;
end
% Elimino appoggio più vecchio di un mese
RawDate3 = find(Dati(:,1)<now-30);
if isempty(RawDate3) == 0
    [rDati,~] = size(Dati);
    if RawDate3(end) == rDati
    else
        Dati = Dati(RawDate3(end)+1:end,:);
    end
end
if isfile(FileTemperature) == 1
    delete(FileTemperature);
end
Dati(:,1) = Dati(:,1) - 730000;
csvwrite(FileTemperature,Dati);

TempDef_IPLHR = TempDef_IPLHR';

fprintf(fileID,fmt,textT);
fclose(fileID);

% azzeramenti di alcuni nodi in particolare
[dsX,dsY,dsZ] = azzeramenti_TLHR(IDcentralina,DTcatena,dsX,dsY,dsZ,NodoInPlaceLinkHR,FileName);

% filtro
[dsX,dsY,dsZ] = filtro(dsX,dsY,dsZ,Ndevst_HR,Wdevst_HR,FileName);
if NuovoZeroIPLHR == 1 && Ndevst_HR ~= 0
    if NdatiMedia > Ndatidespike
        NdatiF = NdatiMedia;
    else
        NdatiF = Ndatidespike;
    end
    ini = round(NdatiF/2)+1;
    if rem(NdatiF,2) == 0
        ini = ini+1;
    end
    clear NDatiF
    iniST = round(Wdevst_HR/2);
    if rem(Wdevst_HR,2) == 0
        iniST = iniST+1;
    end
    iniST = iniST + margine;
    if iniST > ini
        ini = iniST;
    end
    if ini < 6
        ini = 6;
    end
    dsX = dsX(:,ini:end);
    dsY = dsY(:,ini:end);
    dsZ = dsZ(:,ini:end);
    TempDef_IPLHR = TempDef_IPLHR(ini:end,:);
    DatiElabInPlaceLinkHR = DatiElabInPlaceLinkHR(ini:end,:);
    ARRAYdateIPLHR = ARRAYdateIPLHR(ini:end,1); 
    ErrInPlaceLinkHR = ErrInPlaceLinkHR(ini:end,:);
end

%% Finalizzo i calcoli
[rX,cX] = size(dsX);
sommaX_IPLHR = zeros(rIPLHR,1);
Xlocal_IPLHR = zeros(rX,cX+1); % totale del singolo nodo data per data
AlfaX_IPLHR = zeros(rX,cX+1); % Angoli
sommaY_IPLHR = zeros(rIPLHR,1);
Ylocal_IPLHR = zeros(rX,cX+1);
AlfaY_IPLHR = zeros(rX,cX+1); % Angoli
sommaZ_IPLHR = zeros(rX,cX+1);
Zlocal_IPLHR = zeros(rX,cX+1);
X_IPLHR = zeros(rX,cX+1); % cumulata del singolo nodo data per data
Y_IPLHR = zeros(rX,cX+1); % il primo valore è 0
Z_IPLHR = zeros(rX,cX+1);
HShift_IPLHR = zeros(rX,cX+1);
HShift_local_IPLHR = zeros(rX,cX+1);
Azimuth_IPLHR = zeros(rX,cX+1);
azim_IPLHR = zeros(rX,cX+1); % matrice di appoggio per il calcolo dell'azimuth
Speed_IPLHR = zeros(rIPLHR,cX+1); % Velocità 2D Cumulata
Speed_local_IPLHR = zeros(rIPLHR,cX+1); % Velocità 2D locale
Acceleration_IPLHR = zeros(rIPLHR,cX+1); % Accelerazione 2D Cumulata
Acceleration_local_IPLHR = zeros(rIPLHR,cX+1); % Accelerazione 2D Locale
if  yesIPL == 1
    if yesIPLHR3D == 1
        SpeIPLHR = SpeIPLHR(2:end);  % correzione che salta l'ancora
    else
        % Correzione che salta l'ancora e il primo nodo Tilt Link (se
        % presente)
        SpeIPLHR = SpeIPLHR(3:end);
    end
else
    SpeIPLHR = SpeIPLHR(2:end);  % correzione che salta l'ancora
end
% Recupero i dati già elaborati
if NuovoZeroIPLHR == 1
   [rE,cE] = size(DatiElabInPlaceLinkHR);
   cont = 3;
   n = 1;
   while cont<=cE
         sommaX_IPLHR(n,1) = cell2mat(DatiElabInPlaceLinkHR(1,cont))';
         Xlocal_IPLHR(n,1) = sommaX_IPLHR(n,1);
         X_IPLHR(n,1) = cell2mat(DatiElabInPlaceLinkHR(1,cont+3))';
         sommaY_IPLHR(n,1) = cell2mat(DatiElabInPlaceLinkHR(1,cont+1))';
         Ylocal_IPLHR(n,1) = sommaY_IPLHR(n,1);
         Y_IPLHR(n,1) = cell2mat(DatiElabInPlaceLinkHR(1,cont+4))';
         for j = 1:rIPLHR
             AlfaX_IPLHR(j,1) = asind(Xlocal_IPLHR(j,1)/SpeIPLHR(j));
             AlfaY_IPLHR(j,1) = asind(Ylocal_IPLHR(j,1)/SpeIPLHR(j));
         end
         Zlocal_IPLHR(n,1) = cell2mat(DatiElabInPlaceLinkHR(1,cont+2))';
         Z_IPLHR(n,1) = cell2mat(DatiElabInPlaceLinkHR(1,cont+5))';
         HShift_IPLHR(n,1) = cell2mat(DatiElabInPlaceLinkHR(1,cont+6))';
         HShift_local_IPLHR(n,1) = cell2mat(DatiElabInPlaceLinkHR(1,cont+7))';
         Azimuth_IPLHR(n,1) = cell2mat(DatiElabInPlaceLinkHR(1,cont+8))';
         Speed_IPLHR(n,1:rE) = cell2mat(DatiElabInPlaceLinkHR(:,cont+10))';
         Speed_local_IPLHR(n,1:rE) = cell2mat(DatiElabInPlaceLinkHR(:,cont+11))';
         Acceleration_IPLHR(n,1:rE) = cell2mat(DatiElabInPlaceLinkHR(:,cont+12))';
         Acceleration_local_IPLHR(n,1:rE) = cell2mat(DatiElabInPlaceLinkHR(:,cont+13))';
         cont=cont+16;
         n = n+1;
   end
else
    Zlocal_IPLHR(:,1) = SpeIPLHR;
    Z_IPLHR(:,1) = PsIPLHR(2:end);
end

% elaboro i dati nuovi
for iii = 1:cX
    Xlocal_IPLHR(:,iii+1) = sum(dsX(:,1:iii),2)+sommaX_IPLHR(:,1);
    Ylocal_IPLHR(:,iii+1) = sum(dsY(:,1:iii),2)+sommaY_IPLHR(:,1);
    for j = 1:rIPLHR
        AlfaX_IPLHR(j,iii+1) = asind(Xlocal_IPLHR(j,iii+1)/SpeIPLHR(j));
        AlfaY_IPLHR(j,iii+1) = asind(Ylocal_IPLHR(j,iii+1)/SpeIPLHR(j));
    end
    sommaZ_IPLHR(:,iii+1) = sum(dsZ(:,1:iii),2);
    X_IPLHR(:,iii+1) = cumsum(Xlocal_IPLHR(:,iii+1));
    Y_IPLHR(:,iii+1) = cumsum(Ylocal_IPLHR(:,iii+1));
    Z_IPLHR(:,iii+1) = Z_IPLHR(:,1) - cumsum(sommaZ_IPLHR(:,iii+1));
    HShift_IPLHR(:,iii+1) = (X_IPLHR(:,iii+1).^2+Y_IPLHR(:,iii+1).^2).^0.5;
    HShift_local_IPLHR(:,iii+1) = (Xlocal_IPLHR(:,iii+1).^2+Ylocal_IPLHR(:,iii+1).^2).^0.5;
    Zlocal_IPLHR(:,iii+1) = SpeIPLHR - sommaZ_IPLHR(:,iii+1); % Zeta è il singolo abbassamento di quel nodo
    for rr = 1:rIPLHR
        azim_IPLHR(rr,iii) = (acos(abs(X_IPLHR(rr,iii))/HShift_IPLHR(rr,iii)))*180/3.141592654; % Angolo Teta in gradi 0° - 90°
        segnoNS = sign(X_IPLHR(rr,iii));
        segnoEO = sign(Y_IPLHR(rr,iii));
        % L'azimuth si calcola con NS = 0, positivo in senso orario
        % (90° = Est)
        if segnoNS == 1 && segnoEO == 1 % quadrante 1
            Azimuth_IPLHR(rr,iii+1) = azim_IPLHR(rr,iii); % Teta lo tengo come è (1 quadrante)
        elseif segnoNS == -1 && segnoEO == 1 % quadrante 2
               Azimuth_IPLHR(rr,iii+1) = 180 - azim_IPLHR(rr,iii); % 180-teta
        elseif segnoNS == -1 && segnoEO == -1 % quadrante 3
               Azimuth_IPLHR(rr,iii+1) = 180 + azim_IPLHR(rr,iii); % 180+teta
        elseif segnoNS == 1 && segnoEO == -1 % quadrante 4
               Azimuth_IPLHR(rr,iii+1) = 360 - azim_IPLHR(rr,iii); % 360-teta
        elseif segnoNS == 0 && segnoEO == -1 % 270°
               Azimuth_IPLHR(rr,iii+1) = 270;
        elseif segnoNS == -1 && segnoEO == 0 % 180°
               Azimuth_IPLHR(rr,iii+1) = 180;
        end
    end
end
Azimuth_IPLHR = real(Azimuth_IPLHR);
Azimuth_IPLHR = Azimuth_IPLHR+Corr_Azimuth;
[rAz,cAz] = size(Azimuth_IPLHR);
for a = 1:rAz
    for b = 1:cAz
        if Azimuth_IPLHR(a,b) >= 360
           Azimuth_IPLHR(a,b) = Azimuth_IPLHR(a,b)-360;
        end
    end
end

%% Calcolo velocità di spostamento giornaliera
[numDate,~] = size(ARRAYdateIPLHR); % numero di date
d = 1;
p = 1;
diffDate = 0;
n = 1;
period = 1; % calcolo giornaliero
ContSUP = [];
for dd = 1:numDate
    while diffDate < period
          d = d+1;  
          if d > numDate % Se d supera le date disponibili, esco dal ciclo while
              break
          end
          diffDate = ARRAYdateIPLHR(d) - ARRAYdateIPLHR(p);
    end
    if d >numDate
        break
    end
    ContSUP(n,1) = d; % Creo matrice indici dell'estremo superiore della differenza
    ContINF(n,1) = p; % Creo matrice indici dell'estremo inferiore della differenza
    p = p+1; % passo alla data di partenza successiva 
    d = p; % resetto il conto di d
    n = n+1;
    diffDate = 0;
end

check = isempty(ContSUP);
if check == 0
    [nDate,~] = size(ContSUP);
    for s = 1:rIPLHR
        N = 1;
        for dd = 1:nDate
            Speed_IPLHR(s,ContSUP(N,1)) = (HShift_IPLHR(s,ContSUP(N,1))-HShift_IPLHR(s,ContINF(N,1)))/(ARRAYdateIPLHR(ContSUP(N,1))-ARRAYdateIPLHR(ContINF(N,1)));
            Speed_local_IPLHR(s,ContSUP(N,1)) = (HShift_local_IPLHR(s,ContSUP(N,1))-HShift_local_IPLHR(s,ContINF(N,1)))/(ARRAYdateIPLHR(ContSUP(N,1))-ARRAYdateIPLHR(ContINF(N,1)));
            Acceleration_IPLHR(s,ContSUP(N,1)) = (Speed_IPLHR(s,ContSUP(N,1))-Speed_IPLHR(s,ContINF(N,1)))/(ARRAYdateIPLHR(ContSUP(N,1))-ARRAYdateIPLHR(ContINF(N,1)));
            Acceleration_local_IPLHR(s,ContSUP(N,1)) = (Speed_local_IPLHR(s,ContSUP(N,1))-Speed_local_IPLHR(s,ContINF(N,1)))/(ARRAYdateIPLHR(ContSUP(N,1))-ARRAYdateIPLHR(ContINF(N,1)));
            N = N+1;
        end
    end
end

% Approssimo i dati con il corretto numero di cifre decimali
[X_IPLHR,Y_IPLHR,Z_IPLHR,Xlocal_IPLHR,Ylocal_IPLHR,Zlocal_IPLHR,HShift_IPLHR,...
    HShift_local_IPLHR,Azimuth_IPLHR,Speed_IPLHR,Speed_local_IPLHR,Acceleration_IPLHR,...
    Acceleration_local_IPLHR,TempDef_IPLHR] = approx_TLHR(X_IPLHR,Y_IPLHR,Z_IPLHR,...
    Xlocal_IPLHR,Ylocal_IPLHR,Zlocal_IPLHR,HShift_IPLHR,HShift_local_IPLHR,Azimuth_IPLHR,...
    Speed_IPLHR,Speed_local_IPLHR,Acceleration_IPLHR,Acceleration_local_IPLHR,TempDef_IPLHR,FileName);

% Matrice errori
[r,~] = size(ErrInPlaceLinkHR);
Matrice_err = zeros(r,rIPLHR);
if strcmp(NodoInPlaceLinkHR(1,4),'ADC') || strcmp(NodoInPlaceLinkHR{1,4},'null') || isempty(NodoInPlaceLinkHR{1,4}) == 1
    col = 6;
elseif strcmp(NodoInPlaceLinkHR{1,4} ,'digit')
    col = 3;
end
for i = 1:r % date
    d = 1;
    for n = 1:rIPLHR % nodi
        j = 1;
        err = ErrInPlaceLinkHR(i,d:d+col-1);
        while j <= col
            if err(1,j) == 1
                Matrice_err(i,n) = 1;
                break
            end
            if err(1,j) == 0.5
                Matrice_err(i,n) = 0.5;
            end
            j = j+1;
        end
        d = d+col;
    end
end
ErrInPlaceLinkHR = Matrice_err';

text = 'In Place Link HR calculation executed correctly';
fileID = fopen(FileName,'a');
fmt = '%s \r';
fprintf(fileID,fmt,text);
fclose(fileID);

end