view src/mz_hdf5tosci.sci @ 0:0d2d20ef6e72 draft

Uploaded
author jcb-mpl
date Tue, 27 Apr 2021 13:52:21 +0000
parents
children
line wrap: on
line source

function x_mz=mz_hdf5tosci(file_in,max_scans)
    
    // entrees: -------------------------
    // file_in: un fichier au format HDF5 
    // max_scans: valeur introduite pour accelerer la fonction de test
    
    // sorties: ---------------------------
    // x_mz: une variable avec deux champs:
    // x_mz.time
    // x_mz.mzdata une liste dont chaque élément est un vecteur de 2 colonnes: m/z et intensités de signal
    
    // JCB 2 aout18
    
    h5close()
    
    // lecture du fichier HDF5 ----------------------------------
    a=h5open(file_in,'r','stdio');  // que lecture   6jan20:stdio 
    //pause
    spectrum_chromatogram_index=h5read(a.root.ChromatogramIndex);
    spectrum_chromatogram_index=spectrum_chromatogram_index(1);    // une valeur 
    spectrum_time=h5read(a.root.ChomatogramTime);
    spectrum_intensity=h5read(a.root.SpectrumIntensity);
    spectrum_mz=h5read(a.root.SpectrumMZ);
    spectrum_index=h5read(a.root.SpectrumIndex);                    
    //spectrum_index=spectrum_index(1:spectrum_chromatogram_index); 
    h5close(a)
    
    // reconstruction du spectre -----------------
    if argn(2)<2 then
        nbr_scans=max(size(spectrum_index));
    else
        nbr_scans=max_scans;
    end

 
    // 1- extraction du temps
    //25juin20
    time=spectrum_time;   
    if size(time,1)==1 then
        time=time';
    end
    //pas=(spectrum_time($)-spectrum_time(1))/(nbr_scans-1);
    //time=[spectrum_time(1):pas:spectrum_time($)+pas]; // rajout de 1 par securite pour être sur d'arriver a nbr_scans 
    //time=time(1:nbr_scans);    // nbr exact de RT 

    mzdata=list();
    
    //disp(nbr_scans,'nbr_scans=')
    
    for i=1:nbr_scans;
        if i==1 then 
            plage_i=[1,double(spectrum_index(i))];
        else 
            plage_i=[double(spectrum_index(i-1))+1,double(spectrum_index(i))];
        end
        //disp(plage_i,'plage_i=')
        // 2- extraction et reconstruction des m/z
        // dans le fichier HDF5, ne sont enregistrées que les différences de m/z entre 2 mesures
        mz_i=spectrum_mz(plage_i(1):plage_i(2));
        
        n=size(mz_i,2);
        mz_i2=zeros(1,n);
        
        mz_i2(1)=mz_i(1);
        for j=2:n;
            mz_i2(j)=mz_i2(j-1)+mz_i(j);
        end
        
        // 3- extraction des signaux
        sig_i=spectrum_intensity(plage_i(1):plage_i(2));

        // regroupement des m/z et signaux
        sp_i=[mz_i2;sig_i];
        
        mzdata(i)=sp_i';
       
        
    end
    
    // verif. que time est un vecteur colonne
    if size(time,2)>1 then
        time=time';
    end
    

    // amelioration de la coherence des donnees 
    index=time(2:$)-time(1:$-1);
    index2=find(index<0);
    if index2 ~=[] then
        index3=min(index2); 
    else;
        index3=max(size(time));
    end
    
    x_mz.time=time(1:index3);
    mzdata2=list();
    for i=1:min(index3,size(mzdata));             
        mzdata2(i)=mzdata(i);
    end
    x_mz.mzdata=mzdata2;
    
    
endfunction