
### The script is very long, but most parts are skipped with if conditions

 pkg load  spline-gcvspl
 c_colori;  ## a script which contains personal color definitions

pathw="/home/lidia/Tempwork/";
pathd='/home/lidia/Tempwork/';

cd(pathw);

ff_v=5; fig_v=[12]; xi=-24; xf=36;


gfmese_NB=[0,31,59,90,120,151,181,212,243,273,304,334,365];
gfmese_BB=[0,31,60,91,121,152,182,213,244,274,305,335,366];


if(!exist('filess'))
  stringa='*c02';
  if(!exist('stringa'))
    stringa='';
  endif
  filend='*info.dat';

  cd(pathd)
  [stato,filesls]=system(['ls  '  stringa filend]);
  cd(pathw);
  clear aa; aa=strfind(filesls,"\n");
  filess=filesls(1:aa(1)-1);
  for ii=2:length(aa)
   filess=[filess; filesls(aa(ii-1)+1:aa(ii)-1)];
  endfor
endif


evsi=1;
evento_sid='2011 Giappone';

if(evsi==1 & !exist('evento_sid'))
   'File con maremoto? specifica evento_sid. Se no, evsi!=1. STOP.'
   Stop
else
   evento_sid='2011 Giappone'                      ;
endif


"Attenzione valori da det o cambiare: "                   ;
filterfirsi=0                   ;
calcdetezsi=0                   ;
figuresi=1                   ;
fig2printsi=0;
fig4printsi=0;

'da det: fig_v e caso_v'                   ;
'file con la data: yy mm dd hh mm ss'                    ;
ev_te_ss=0;


if(!exist("printsi"))
  printsi=0                   ;
endif
if(!exist("filterfirsi"))
    "filter fir ? filterfirsi=1 (sì). Det filterfirsi. Stop."
    filterfirsi=0;
endif
if(!exist("calcdetezsi"))
    "calcola detezioni ? calcdetezsi=1 (sì). Det calcdetezsi. Stop."
    calcdetezsi=0                   ;
endif
if(!exist("calcslopeteorsi"))
   calcslopeteorsi=0                   ;
endif
if(!exist("ff_v") & exist('ff'))
  ff_v=ff;
endif
if(!exist("checksi"))
  checksi=1456;
endif
if(!exist("new_smoothsi"))
  new_smoothsi=-56;
endif
if(!exist("ffp"))
  ffp=0;
endif
if(!exist("dett"))
  dett=[];
endif
if(!exist("ccp"))
  ccp=-1;
endif
if(!exist("solofig2si"))
  solofig2si=-1;
endif
if(!exist("fig8si"))
  fig8si=-1;
endif
if(!exist('fig_v'))
  fig_v=[];
endif
if(fig8si==1 & calcslopeteorsi==1)
   fig_v=[fig_v; 8];
endif
if(figuresi==1)
  for ii=fig_v
     figure(ii); clf;
  endfor
endif


nf=size(filess,1);
nMcasi=7;
nf_ta=[];  nf_calc=[];  nf_calc_sm=[];  nt_ta=[];  nt_calc=[];  nt_calc_sm=[];
stima_mareaa_MM=[]; stima_mareaa_mm=[]; stima_mareaa_md=[];


caso_v=2;



ris=[];
c_iis=[             ### personal colors definition from c_colori
       c_verde;
       c_verde;
       c_azzurro1;
       c_azzurro1;
       c_blu;
       c_blu;
       c_viola3;];



for ff=ff_v
for caso=caso_v

  if(!exist("scalatemp"))
     scalatemp="hh";
     sst=1/60;
  elseif(scalatemp=="gg")
     sst=1/1440;
  elseif(scalatemp=="hh")
     sst=1/60;
  elseif(scalatemp=="mm")
     sst=1;
  elseif(scalatemp=="ss")
     sst=60;
  else
     sst=1;
  endif

  if(ffp!=ff | !exist("y_sl") | ccp<=0 | ccp!=caso)
'loading'
    clear -x c_* ccp dett ff* caso* checksi fig2printsi xi xf gfmese* path* filess printsi scalatemp nMcasi sst calcdetezsi   figuresi nf_* nt_* stima_mareaa* T* filterfirsi fig8si solofig2si N kk kk_v N_v ris fig_v ll* leg* new_smoothsi calcslopeteorsi ev_te_ss evsi evento_sid stringa fig2printsi fig4printsi tat ENM k_filtro ert asi;

    if(caso<10)
      casos=["0" num2str(caso)];
    else
      casos=num2str(caso);
    endif

    tfile=filess(ff,:);
    radfile=tfile(1:index(tfile,"_c")-1)                   ;
    sito=tfile(1:index(tfile,"_")-1);

    file_TF = [ radfile  "_c" casos "_costrIS.dat"];
    fileTD  = [ radfile  "_c" casos "_ISBSCF.dat"]                   ;
    fileMsd = [ radfile  "_c" casos "_MSDvar.dat"];
    fileInf = [ radfile  "_c" casos "_tainfo.dat"];
    source([ pathd fileInf]);

    dat_TF = load([ pathd file_TF ]);
    datTD  = load([ pathd  fileTD  ]);
    datMsd = load([ pathd  fileMsd ]);

    dat=dat_TF;
    ndat=size(dat,1);

    anno=dat(:,1);
    iud=6;
    if(max(diff(anno))==0)
       clear bb; bb=find(mod(anno,4)==0  & (!(mod(anno,100)==0 & mod(anno,400)!=0)));
       clear nb; nb=find(!(mod(anno,4)==0  & (!(mod(anno,100)==0 & mod(anno,400)!=0))));
       if(!isempty(nb))
          xxd(nb)=((gfmese_NB(dat(nb,2))'.+dat(nb,3)).*24.+dat(nb,4)).*60.+dat(nb,5).+dat(nb,6)./60;
       endif
       if(!isempty(bb))
          xxd(bb)=((gfmese_BB(dat(bb,2))'.+dat(bb,3)).*24.+dat(bb,4)).*60.+dat(bb,5).+dat(bb,6)./60;
       endif
    else
       gg_annoprec=0;
       for annox=min(anno):max(anno)
           clear yy; yy=find(dat(:,1)==annox);
           if(mod(annox,4)==0  & (!(mod(annox,100)==0 & mod(annox,400)!=0)))
              xxd(yy)=((gg_annoprec.+gfmese_BB(dat(yy,2))'.+dat(yy,3)).*24.+dat(yy,4)).*60.+dat(yy,5).+dat(yy,6)./60;
              gg_annoprec=gg_annoprec+max(gfmese_BB);
           else
              xxd(yy)=((gg_annoprec.+gfmese_NB(dat(yy,2))'.+dat(yy,3)).*24.+dat(yy,4)).*60.+dat(yy,5).+dat(yy,6)./60;
              gg_annoprec=gg_annoprec+max(gfmese_NB);
           endif
       endfor
    endif

    if(ev_te_ss==1)
       data_rif=[0 1 0 0 0 0 ];
    elseif(evsi!=1)
       data_rif=[dat(1,1:2) 1 0 0 0 ];
    elseif(evsi==1)
       dati_eventi_source;  ## script with date information
       for ii=1:size(origine,1)
         if(index(origine,evento_sid)==1)
           ind=ii;
         else
           continue
         endif
       endfor
       date_eq=data_eq(ind,:);  #### date_eq = [2011      3     11      5     46     23]
       data_rif=date_eq;
    endif
    anno=data_rif(1);
    if(mod(min(anno),4)==0  & (!(mod(min(anno),100)==0 & mod(min(anno),400)!=0)))
       pshift=((gfmese_BB(data_rif(2))'.+data_rif(3)).*24.+data_rif(4)).*60.+data_rif(5).+data_rif(6)./60-1;
    else
       pshift=((gfmese_NB(data_rif(2))'.+data_rif(3)).*24.+data_rif(4)).*60.+data_rif(5).+data_rif(6)./60-1;
    endif
    xxd=xxd'.-pshift;

    if(checksi==1)
        clear dat xxt; dat=datTD;
        xxt=(((gfmese(dat(:,2))'.+dat(:,3)).*24).+dat(:,4)).*60.+dat(:,5).+dat(:,6)./60;  ## xxd in minuti
        xxt=xxt.-pshift;
        if(max(abs(xxd.-xxt))>10^-6)
            ErroreeeeeeSTOP
        endif
        clear dat xxt; dat=datMsd;
        xxt=(((gfmese(dat(:,2))'.+dat(:,3)).*24).+dat(:,4)).*60.+dat(:,5).+dat(:,6)./60;  ## xxd in minuti
        xxt=xxt.-xxt(1);
        if(max(abs(xxd.-xxt))>10^-6)
            ErroreeeeeeSTOP
        endif
    endif

    y_slt=dat(:,iud+1);
    IST  =dat(:,iud+2);
    tf_uns =dat(:,iud+3);
    TF   =dat(:,iud+4);
    IS   =dat(:,iud+5);
    clear dat dat_TF;

    dat=datTD;
    y_sl =dat(:,iud+1);
    f_IS =dat(:,iud+2);
    f_BS1=dat(:,iud+3);
    f_CF1=dat(:,iud+4);
    f_BS2=dat(:,iud+5);
    f_CF2=dat(:,iud+6);
    f_BS3=dat(:,iud+7);
    f_CF3=dat(:,iud+8);
    if(max(abs(y_sl.-y_slt))>10^-8)
        ErroreeeeeeSTOP
    else
       clear y_slt
    endif
    if(max(abs(IS.-f_IS))>10^-8)
        ErroreeeeeeSTOP
    else
        clear IS
    endif
    clear dat datTD2;

    dat=datMsd;
    y_slt=dat(:,iud+1);
    Msd  =dat(:,iud+2);
    varz =dat(:,iud+3);
    v_sp =dat(:,iud+4);
    v_ns =dat(:,iud+5);
    varR =dat(:,iud+6);
    if(max(abs(y_sl.-y_slt))>10^-8)
        ErroreeeeeeSTOP
    else
     clear y_slt
    endif
    clear dat datMsd;

    if(!exist('k_filtro'))
      k_filtro=0.000001;
    endif
    tide_fy_sl    = csaps(xxd, y_sl,k_filtro, xxd);
    y_slfil = y_sl .- tide_fy_sl;


# # # # filtro a posteriori
# # # # -------------------
    if(calcslopeteorsi==1)
       sl_y_sl_teor=[0; diff(y_sl)];
       sl_tide_fy_sl_teor=[0; diff(tide_fy_sl)];
       sl_y_slfil_teor=[0; diff(y_slfil)];
       if(abs(sl_y_slfil_teor.-(sl_y_sl_teor.-sl_tide_fy_sl_teor))>10^-8)
           ErroreeeeeeSTOP
       endif
       for ii=Nt_IS:ndat
         y_slfil_sm(ii)=mean(y_slfil(ii-Nt_IS+1:ii));
         sl_y_sl_teor_sm(ii)=mean(sl_y_sl_teor(ii-Nt_IS+1:ii));
         sl_y_slfil_teor_sm(ii)=mean(sl_y_slfil_teor(ii-Nt_IS+1:ii));
       endfor
    endif

    casop=caso;
    ffp=ff;
    ccp=caso;

'loading done'
  endif

# # # # -----------------
  if(new_smoothsi==1)
     Nt_Smooth=floor(t_Smooth_f*60/samplingrate)
     if(Nt_Smooth>0)
       clear tf2;
       tmm=max(60,Nt_Smooth);
       for ii=1:Nt_Smooth-1
         tf2(ii)=tf_uns(ii);
       endfor
       for ii=Nt_Smooth:tmm-1
         tf2(ii)=tf_uns(ii);
         tf2(ii)=mean(tf2(ii-Nt_Smooth+1:ii));
       endfor
       for ii=tmm:ndat
         tf2(ii)=tf_uns(ii);
         tf2(ii)=mean(tf2(ii-Nt_Smooth+1:ii));
       endfor
       tf2=tf2';
     else
       tf2=TF;
     endif
  else
     tf2=TF;
  endif
# # # # -----------------

  if(!exist("xi"))
    xi=xxd(1).*sst;
  endif
  if(!exist("xf"))
    xf=xxd(end)*sst;
  endif




  iicd=Nt_IS.+Nt_TF.+Nt_G_TF.+Nt_BS.+Nt_G_BS.-2;



  if(filterfirsi==1)
# # # #        [n, Wn, bet, ftype] = kaiserord(f, m, dev [, fs]);
# # # #        n: order of the filter (length of filter minus 1)
# # # #        Wn: band edges for use in fir1
# # # #        beta: parameter for kaiser window of length n+1
# # # #        ftype: choose between pass and stop bands
# # # #        bfir = fir1(n,Wn,kaiser(n+1,beta),ftype,'noscale');
# # # #         N=16; bfir  = fir1(N,2/(3*60) ,"stop","noscale" ); freqz(bfir);
# # # #         bfir  = fir1(N,1/(3*60) ,"stop","noscale" );
# # # #         bfir  = fir1(N, 2/(3*60) ,"stop","noscale" );
       filtr=1;
        Nfir=200;
        bfir  = fir1(Nfir,2/(3*60) ,"high","kaiser","noscale" ); ## freqz(bfir);
        bfir  = fir1(Nfir,2/(3*60) ,"high","noscale" ); ## freqz(bfir);
        y_fir=y_sl.*0;
        y_fir(1:Nfir+1)=filter(bfir,1,y_sl(1:Nfir+1));
        y_fir(1:Nfir)=y_fir(Nfir+1).*[1:Nfir]./(Nfir+1);
        for ii=Nfir+2:ndat
           y_fir_tot = filter(bfir,1,y_sl(ii-Nfir:ii));
           y_fir(ii)=y_fir_tot(Nfir+1);
        endfor
       sl_fir=[0; diff(y_fir)];
       sl_fir_sm=sl_fir.*0;
       for ii=Nt_IS:ndat
         sl_fir_sm(ii)=mean(sl_fir(ii-Nt_IS+1:ii));
       endfor

       ris=[ ris; [Nfir kk filtr abs(max(y_slfil)-max(y_fir))   abs(min(y_slfil)-min(y_fir))  abs(((max(y_slfil)-min(y_slfil))-(max(y_fir)-min(y_fir)))) ]];

# # # #  Butterworth filter
# # # #  [b,a] = butter(n, Wc, 'high')
# # # #  high pass filter with cutoff pi*Wc radians
# # # # --------------------------------------------
# # # # N_v=[1:13];
# # # # kk_v=[1 1/pi];
       filtr=2;
       kk=1/pi; Nbutt=14;
       [bbutt,abutt] = butter(Nbutt, kk*2/(3*60), 'high');
       y_butt=y_sl.*0;
       y_butt(1:Nbutt+1)=filter(bbutt,abutt,y_sl(1:Nbutt+1));
       y_butt(1:Nbutt)=y_butt(Nbutt+1).*[1:Nbutt]./(Nbutt+1);
       for ii=Nbutt+2:ndat
          y_butt_tot = filter(bbutt,abutt,y_sl(ii-Nbutt:ii));
          y_butt(ii)=y_butt_tot(Nbutt+1);
       endfor
       sl_butt=[0; diff(y_butt)];
       sl_butt_sm=sl_butt.*0;
       for ii=Nt_IS:ndat
         sl_butt_sm(ii)=mean(sl_butt(ii-Nt_IS+1:ii));
       endfor

       ris=[ ris; [Nbutt kk filtr abs(max(y_slfil)-max(y_butt))   abs(min(y_slfil)-min(y_butt))  abs(((max(y_slfil)-min(y_slfil))-(max(y_butt)-min(y_butt)))) ]];

  endif



  if(calcdetezsi==1)
    for met=1:3
       if(met==3 & caso==5)
          clear aa; aa=find(RangeDDET(:,1)==met & RangeDDET(:,2)==caso & RangeDDET(:,5)==5);
       else
          clear aa; aa=find(RangeDDET(:,1)==met & RangeDDET(:,2)==caso);
       endif
       if(!(RangeDDET(aa,3)==5 & RangeDDET(aa,4)==-5))
          lCF_v(met,:)=[RangeDDET(aa,3) RangeDDET(aa,4)];
       else
          lCF_v(met,:)=[5 5];
       endif
       if(max(size(aa))>1)
          ErroreeeeeeSTOP
       endif
    endfor
    clear det_v1m_ta;    det_v1m_ta   =find(abs(IS)>=1 & abs(IS)>=lCF_v(1,1).*BS1 & xxd>=xxd(iicd));
    clear det_v1M_ta;    det_v1M_ta   =find(abs(IS)>=1 & abs(IS)>=lCF_v(1,2).*BS1 & xxd>=xxd(iicd));
    clear det_v2m_ta;    det_v2m_ta   =find(abs(IS)>=1 & abs(IS)>=lCF_v(2,1).*BS2 & xxd>=xxd(iicd));
    clear det_v2M_ta;    det_v2M_ta   =find(abs(IS)>=1 & abs(IS)>=lCF_v(2,2).*BS2 & xxd>=xxd(iicd));
    clear det_v3m_ta;    det_v3m_ta   =find(abs(IS)>=1 & abs(IS)>=lCF_v(3,1).*BS3 & xxd>=xxd(iicd));
    clear det_v3M_ta;    det_v3M_ta   =find(abs(IS)>=1 & abs(IS)>=lCF_v(3,2).*BS3 & xxd>=xxd(iicd));

    nt_v1m_ta=size(det_v1m_ta,1);    nF_v1m_ta=size(find(xxd(det_v1m_ta)<0 | xxd(det_v1m_ta)>durtsu),1);
    nt_v1M_ta=size(det_v1M_ta,1);    nF_v1M_ta=size(find(xxd(det_v1M_ta)<0 | xxd(det_v1M_ta)>durtsu),1);
    nt_v2m_ta=size(det_v2m_ta,1);    nF_v2m_ta=size(find(xxd(det_v2m_ta)<0 | xxd(det_v2m_ta)>durtsu),1);
    nt_v2M_ta=size(det_v2M_ta,1);    nF_v2M_ta=size(find(xxd(det_v2M_ta)<0 | xxd(det_v2M_ta)>durtsu),1);
    nt_v3m_ta=size(det_v3m_ta,1);    nF_v3m_ta=size(find(xxd(det_v3m_ta)<0 | xxd(det_v3m_ta)>durtsu),1);
    nt_v3M_ta=size(det_v3M_ta,1);    nF_v3M_ta=size(find(xxd(det_v3M_ta)<0 | xxd(det_v3M_ta)>durtsu),1);

    iidetm_v1m_ta   =find(xxd(det_v1m_ta   )>=0 & xxd(det_v1m_ta)<=durtsu);
    iidetm_v1M_ta   =find(xxd(det_v1M_ta   )>=0 & xxd(det_v1M_ta)<=durtsu);
    iidetm_v2m_ta   =find(xxd(det_v2m_ta   )>=0 & xxd(det_v2m_ta)<=durtsu);
    iidetm_v2M_ta   =find(xxd(det_v2M_ta   )>=0 & xxd(det_v2M_ta)<=durtsu);
    iidetm_v3m_ta   =find(xxd(det_v3m_ta   )>=0 & xxd(det_v3m_ta)<=durtsu);
    iidetm_v3M_ta   =find(xxd(det_v3M_ta   )>=0 & xxd(det_v3M_ta)<=durtsu);

    nf_ta   =[ nf_ta   ; [ff  nF_v1m_ta    nF_v1M_ta    nF_v2m_ta    nF_v2M_ta    nF_v3m_ta    nF_v3M_ta ]    ];
    nt_ta   =[ nt_ta   ; [ff  nt_v1m_ta    nt_v1M_ta    nt_v2m_ta    nt_v2M_ta    nt_v3m_ta    nt_v3M_ta ]   ];

    xdet_v1m_ta   =repmat(xxd(det_v1m_ta   ).*sst,1,2)';
    xdet_v1M_ta   =repmat(xxd(det_v1M_ta   ).*sst,1,2)';
    xdet_v2m_ta   =repmat(xxd(det_v2m_ta   ).*sst,1,2)';
    xdet_v2M_ta   =repmat(xxd(det_v2M_ta   ).*sst,1,2)';
    xdet_v3m_ta   =repmat(xxd(det_v3m_ta   ).*sst,1,2)';
    xdet_v3M_ta   =repmat(xxd(det_v3M_ta   ).*sst,1,2)';

  endif

  if(calcslopeteorsi==1)
    stima_mareaa_MM=[stima_mareaa_MM; [ff   max(abs(tf2.-sl_tide_fy_sl_teor)) max(abs(tf2(mar).-sl_tide_fy_sl_teor(mar))) max(abs(tf2(nomar).-sl_tide_fy_sl_teor(nomar))) ]];
    stima_mareaa_mm=[stima_mareaa_mm; [ff    min(abs(tf2.-sl_tide_fy_sl_teor)) min(abs(tf2(mar).-sl_tide_fy_sl_teor(mar))) min(abs(tf2(nomar).-sl_tide_fy_sl_teor(nomar))) ]];
    stima_mareaa_md=[stima_mareaa_md; [ff    mean(abs(tf2.-sl_tide_fy_sl_teor)) mean(abs(tf2(mar).-sl_tide_fy_sl_teor(mar))) mean(abs(tf2(nomar).-sl_tide_fy_sl_teor(nomar))) ]];
  endif




  if(figuresi==1)

    if(!exist("xi"))
       xi=xxd(1).*sst;
    endif
    if(!exist("xf"))
       xf=xxd(end).*sst;
    endif
    if(xi<xxd(1).*sst-3)
       xi=xxd(1).*sst;
    endif
    if(xf>xxd(end).*sst+3)
       xf=xxd(end).*sst;
    endif

    ijj=find(xxd.*sst>=xi & xxd.*sst<=xf);

    if(!isempty(find(fig_v==1)))
      figure(1);
       hold on; grid on
       plot(xxd(ijj).*sst,y_sl(ijj)      ,";m(t) ;" ,"color",c_blu)
       plot(xxd(ijj).*sst,tide_fy_sl(ijj),";m_d(t);","color",c_verde)
       axis([xi xf])
    endif

    if(!isempty(find(fig_v==2)))
      if(caso==caso_v(1))
        ll=0; leg={};
      endif
      figure(2);
       hold on; grid on
       if(caso==caso_v(1))
        plot(xxd(ijj).*sst,y_slfil(ijj),"color",c_nero, "linewidth",4)
        ll=ll+1; leg{ll}="m_d(t)";
       endif
       if(exist("y_fir"))
         plot(xxd(ijj).*sst,y_fir(ijj) ,"color",c_fuxia,"linewidth",2)
         ll=ll+1; leg{ll}="m_{fir }(t)";
       endif
       if(exist("y_butt"))
         plot(xxd(ijj).*sst,y_butt(ijj),"color",c_rosso,"linewidth",2)
         ll=ll+1; leg{ll}="m_{butt}(t)";
       endif
       axis([xi xf])
    endif

    if(!isempty(find(fig_v==12)))
      if(caso==caso_v(1))
        ll1=0; leg1={};
        ll2=0; leg2={};
      endif
      ym=min(y_slfil(ijj))
      yM=max(y_slfil(ijj))
      tat
      ert
'figure'
      figure(12);
       subplot(2,1,1)
       hold on;
       if(exist('tat') & exist('ert'))
          plot((+ert+tat)*sst.*[1 1],[ym yM]);
          plot((-ert+tat)*sst.*[1 1],[ym yM]);
          plot(tat.*[1 1].*sst,[ym yM],"color",c_verde);  ll1=ll1+1; leg1{ll1}="";
       elseif(exist('tat'))
          plot(tat.*[1 1].*sst,[ym yM],"color",c_verde);  ll1=ll1+1; leg1{ll1}="";
       endif
       if(caso==caso_v(1))
        plot([xi xf],[0 0],'-k')
        plot(xxd(ijj).*sst,y_slfil(ijj),"color",c_nero, "linewidth",4)
        ll1=ll1+1; leg1{ll1}="m_d(t)";
       endif
       axis([xi xf])
       grid on
       subplot(2,1,2)
       ym=min(y_sl(ijj));      yM=max(y_sl(ijj));
       hold on;
       if(exist('tat') & exist('ert'))
          plot(( ert+tat).*[1 1].*sst,[ym yM],"color",c_verde);  ll2=ll2+1; leg2{ll2}="";
          plot((-ert+tat).*[1 1].*sst,[ym yM],"color",c_verde);  ll2=ll2+1; leg2{ll2}="";
          plot(tat.*[1 1].*sst,[ym yM],"color",c_verde);  ll2=ll2+1; leg2{ll2}="";
       elseif(exist('tat'))
          plot(tat.*[1 1].*sst,[ym yM],"color",c_verde);  ll2=ll2+1; leg2{ll2}="";
       endif
       plot(xxd(ijj).*sst,tide_fy_sl(ijj),"color",c_rosso); ll2=ll2+1; leg2{ll2}=";m_d(t)";
       plot(xxd(ijj).*sst,y_sl(ijj)     ,"color",c_blu);   ll2=ll2+1; leg2{ll2}="m(t)";
       axis([xi xf])
       grid on
'end figure'
    endif




    if(!isempty(find(fig_v==3)))
      figure(3)
       hold on; grid on
       plot(xxd(ijj).*sst,IST(ijj) ,";IS_T (t)  ;"  ,"color",c_verde,"linewidth",2)
       plot(xxd(ijj).*sst,tf_uns(ijj),";TF_{uns} (t);","color",c_arancio1,"linewidth",3)
       plot(xxd(ijj).*sst,tf2(ijj) ,";TF (t)   ;"   ,"color",c_arancio2,"linewidth",3)
       axis([xi xf])
       axis([xi xf -1 1])
       xlabel(["time (" scalatemp ")"]);
       ylabel("cm/min")
    endif

    if(!isempty(find(fig_v==4)))
      figure(4)
       hold on; grid on
       if(exist("sl_fir_sm"))
         plot(xxd(ijj).*sst,sl_fir(ijj)    ,"; IS( m_{fir}, t) ;",     "color",c_fuxia,"linewidth",2)
         plot(xxd(ijj).*sst,sl_fir_sm(ijj) ,"; IS_{sm}( m_{fir}, t) ;","color",c_fuxia,"linewidth",1)
       endif
       if(exist("sl_butt_sm"))
         plot(xxd(ijj).*sst,sl_butt(ijj)   ,"; IS( m_{butt}, t) ;",    "color",c_rosso,"linewidth",2)
         plot(xxd(ijj).*sst,sl_butt_sm(ijj),"; IS_{sm}( m_{butt}, t);","color",c_rosso,"linewidth",2)
       endif
       plot(xxd(ijj).*sst,f_IS(ijj)        ,"; IS (t)              ;" ,"color",c_verde,"linewidth",2)
       axis([xi xf])
       axis([xi xf -1 1])
    endif

    if(!isempty(find(fig_v==5)))
      figure(5)
       hold on; grid on
       plot(xxd(ijj).*sst,f_IS(ijj),";IS (t)   ;",   "color",c_blu,"linewidth",1)
       plot(xxd(ijj).*sst,IST(ijj), ";IS_T (t)  ;",  "color",c_verde,"linewidth",2)
       plot(xxd(ijj).*sst,tf_uns(ijj),";TF_{uns} (t);","color",c_arancio1,"linewidth",3)
       plot(xxd(ijj).*sst,tf2(ijj), ";TF (t)   ;",   "color",c_arancio2,"linewidth",3)
       axis([xi xf -1.5 1.5])
       grid on
    endif






    if(!isempty(find(fig_v==6)))
      figure(6); clf
       subplot(3,1,1)
       hold on; grid on
       plot(xxd(ijj).*sst,f_IS(ijj)     ,";IS    ;" ,"color",c_verde,"linewidth",3)
       axis([xi xf])

       subplot(3,1,2)
       hold on; grid on
       plot(xxd(ijj).*sst,f_BS1(ijj)     ,";BS1      ;" ,"color",c_rosso   ,"linewidth",2)
       plot(xxd(ijj).*sst,f_BS2(ijj)     ,";BS2      ;" ,"color",c_verde   ,"linewidth",2)
       plot(xxd(ijj).*sst,f_BS3(ijj)     ,";BS3      ;" ,"color",c_blu     ,"linewidth",2)
       axis([xi xf])

       subplot(3,1,3)
       hold on; grid on
       plot(xxd(ijj).*sst,f_CF1(ijj)     ,";CF1      ;" ,"color",c_rosso   ,"linewidth",2)
       plot(xxd(ijj).*sst,f_CF2(ijj)     ,";CF2      ;" ,"color",c_verde   ,"linewidth",2)
       plot(xxd(ijj).*sst,f_CF3(ijj)     ,";CF3      ;" ,"color",c_blu     ,"linewidth",2)
       axis([xi xf])
    endif

    if(!isempty(find(fig_v==7)))
      lIS=1;
      ll=0; leg={};
      figure(7);
       hold on; grid on
       ym=min(f_IS(ijj));      yM=max(f_IS(ijj));
       if(exist('tat'))
          plot(tat.*[1 1].*sst,[ym yM],"color",c_verde);  ll=ll+1; leg{ll}="";
       endif
       plot(xxd(ijj).*sst,f_IS(ijj)           ,"color",c_viola3,"linewidth",2); ll=ll+1; leg{ll}="IS";
       plot(xxd(ijj([1 end])).*sst, lIS.*[1 1],"color",c_nero    ,"linewidth",1); ll=ll+1; leg{ll}="";
       plot(xxd(ijj([1 end])).*sst,-lIS.*[1 1],"color",c_nero    ,"linewidth",1); ll=ll+1; leg{ll}="";
       plot(xxd(ijj([1 end])).*sst, lIS.*[1 1],"color",c_nero    ,"linewidth",1); ll=ll+1; leg{ll}="";
       plot(xxd(ijj([1 end])).*sst,-lIS.*[1 1],"color",c_nero    ,"linewidth",1); ll=ll+1; leg{ll}="";
       if(calcdetezsi==1)
         for lCF=lCF_v(1,:)
          plot(xxd(ijj).*sst, lCF.*f_BS1(ijj) ,"color",c_rosso   ,"linewidth",1); ll=ll+1; leg{ll}="BS1";
          plot(xxd(ijj).*sst,-lCF.*f_BS1(ijj) ,"color",c_rosso   ,"linewidth",1); ll=ll+1; leg{ll}="";
         endfor
         for lCF=lCF_v(2,:)
          plot(xxd(ijj).*sst, lCF.*f_BS2(ijj) ,"color",c_verde   ,"linewidth",1); ll=ll+1; leg{ll}="BS2";
          plot(xxd(ijj).*sst,-lCF.*f_BS2(ijj) ,"color",c_verde   ,"linewidth",1); ll=ll+1; leg{ll}="";
         endfor
         for lCF=lCF_v(3,:)
          plot(xxd(ijj).*sst, lCF.*f_BS3(ijj) ,"color",c_blu     ,"linewidth",1); ll=ll+1; leg{ll}="BS3";
          plot(xxd(ijj).*sst,-lCF.*f_BS3(ijj) ,"color",c_blu     ,"linewidth",1)
         endfor
         ylim=get(gca,"ylim"); l=abs(ylim(2)-ylim(1)); yM=max(ylim)-0.05*l; ym=min(ylim)+0.05*l;
         plot(xdet_v1m_ta,[xdet_v1m_ta(1,:).*0.+ym; xdet_v1m_ta(1,:).*0.+yM], "-1","linewidth",2); ll=ll+1; leg{ll}="";
       endif
       axis([xi xf -5 5])
    endif

  endif

endfor
endfor







if(printsi==1)
  if(max(size(caso_v))==nMcasi)
    casos="a";
  elseif(max(size(caso_v))>1)
    casos="m";
  endif

  cd(pathw);
  dett=[sito '_t' num2str(t_TF_f) '_tS' num2str(t_Smooth_f)]

  if(!isempty(find(fig_v==1)))
    figure(1)
     title("Ampiezza con tide_fy_sl")
     xlabel(["time (" scalatemp ")"]);
     ylabel("cm")
     legend("off")
     nomefileeps=["fig1_amp_m_ff_" num2str(ff) "_c" casos "_d_" dett ".eps"];
     nomefilepng=["fig1_amp_m_ff_" num2str(ff) "_c" casos "_d_" dett ".png"];
     print(nomefilepng,"-tight","-dpng")
   endif

  if(!isempty(find(fig_v==2)))
    figure(2)
     xlabel(["time (" scalatemp(1) ")"],"fontsize",18);
     ylabel("sea level (cm)","fontsize",18)
     if(fig2printsi==1)
       switch (xf-xi)
        case  3
         dx=1;
        case  4
         dx=1;
        case  6
         dx=1;
        case 12
         dx=2;
        case 24
         dx=3;
        otherwise
         dx=floor((xf-xi)/10);
       endswitch
       xtickv=[xi:dx:xf];
       xtickl={};
       for ii=1:size(xtickv,2); xtickl{ii}=num2str(xtickv(ii)-xi); endfor
       set(gca,"xtick",xtickv,"xticklabel",xtickl,"fontsize",18)
       nomefileeps=["fig2_ampsm_ff_" num2str(ff) "_c" casos "_d_" dett "_p.eps"];
       nomefilepng=["fig2_ampsm_ff_" num2str(ff) "_c" casos "_d_" dett "_p.png"];
       legend(leg)
       legend("location","northwest")
     else
       title("Ampiezza senza tide_fy_sl","fontsize",18)
       nomefileeps=["fig2_ampsm_ff_" num2str(ff) "_c" casos "_d_" dett ".eps"];
       nomefilepng=["fig2_ampsm_ff_" num2str(ff) "_c" casos "_d_" dett ".png"];
       legend("off")
     endif
     print(nomefilepng,"-tight","-dpng")
   endif

  if(!isempty(find(fig_v==3)))
    figure(3)
     title("Pendenza con tide_fy_sl e pendenze di tide_fy_sl")
     xlabel(["time (" scalatemp ")"]);
     ylabel("cm/min")
     legend("off")
     nomefileeps=["fig3_slp_m_ff_" num2str(ff) "_c" casos "_d_" dett ".eps"];
     nomefilepng=["fig3_slp_m_ff_" num2str(ff) "_c" casos "_d_" dett ".png"];
     print(nomefilepng,"-tight","-dpng")
   endif

  if(!isempty(find(fig_v==4)))
    figure(4)
     xlabel(["time (" scalatemp(1) ")"],"fontsize",18);
     ylabel("detided slope (cm/min)","fontsize",18)
     if(fig4printsi==1)
       switch (xf-xi)
        case  3
         dx=1;
        case  4
         dx=1;
        case  6
         dx=1;
        case 12
         dx=2;
        case 24
         dx=3;
        otherwise
         dx=floor((xf-xi)/10);
       endswitch
       xtickv=[xi:dx:xf];
       xtickl={};
       for ii=1:size(xtickv,2); xtickl{ii}=num2str(xtickv(ii)-xi); endfor
       set(gca,"xtick",xtickv,"xticklabel",xtickl,"fontsize",18)
       nomefileeps=["fig4_slpsm_ff_" num2str(ff) "_c" casos "_d_" dett "_p.eps"];
       nomefilepng=["fig4_slpsm_ff_" num2str(ff) "_c" casos "_d_" dett "_p.png"];
     else
       title("Pendenza senza tide_fy_sl")
       legend("off")
       nomefileeps=["fig4_slpsm_ff_" num2str(ff) "_c" casos "_d_" dett ".eps"];
       nomefilepng=["fig4_slpsm_ff_" num2str(ff) "_c" casos "_d_" dett ".png"];
     endif
     print(nomefilepng,"-tight","-dpng")
   endif

  if(!isempty(find(fig_v==5)))
    figure(5)
     switch (xf-xi)
      case  3
       dx=1;
      case  4
       dx=1;
      case  6
       dx=1;
      case 12
       dx=2;
      case 24
       dx=3;
      otherwise
       dx=floor((xf-xi)/10);
     endswitch
     xtickv=[xi:dx:xf];
     xtickl={};
     for ii=1:size(xtickv,2); xtickl{ii}=num2str(xtickv(ii)); endfor
     xlabel(["time (" scalatemp(1) ")"],"fontsize",18);
     ylabel("detided slope (cm/min)","fontsize",18)
     set(gca,"xtick",xtickv,"xticklabel",xtickl,"fontsize",18)
     nomefileeps=["fig5_d_" dett "_p.eps"];
     nomefilepng=["fig5_d_" dett "_p.png"];
     print(nomefilepng,"-tight","-dpng")
   endif



  if(!isempty(find(fig_v==6)))
    figure(6)
     title("Funzioni IS BS e CF")
     xlabel(["time (" scalatemp ")"]);
     ylabel("cm/min")
     legend("off")
     nomefileeps=["fig6_fTEDA_ff_" num2str(ff) "_c" casos "_d_" dett ".eps"];
     nomefilepng=["fig6_fTEDA_ff_" num2str(ff) "_c" casos "_d_" dett ".png"];
     print(nomefilepng,"-tight","-dpng")
   endif


  cd(pathw);
endif




if(fig8si==1 & max(size(ff_v))>1)
  figure(8); clf;
   hold on
   plot(stima_mareaa_MM(:,4),"-1")
   plot(stima_mareaa_MM(:,5),"-3")
   plot(stima_mareaa_md(:,4),"-1")
   plot(stima_mareaa_md(:,5),"-3")
   grid on

endif


if(exist('tat') & exist('ert') & exist('ENM') )
  [radfile '     '  num2str(tat)  '     '  num2str(ert) '               '  num2str(ENM)]
endif


