Changeset 939 for trunk/src


Ignore:
Timestamp:
Apr 6, 2016, 10:35:30 PM (9 years ago)
Author:
sommeria
Message:

extract corrected and merge_proj_polar added

Location:
trunk/src
Files:
4 edited

Legend:

Unmodified
Added
Removed

  • trunk/src/series/extract_multitif.m

    r931 r939  
    6565    ParamOut.AllowInputSort='off';% allow alphabetic sorting of the list of input file SubDir (options 'off'/'on', 'off' by default)
    6666    ParamOut.WholeIndexRange='on';% prescribes the file index ranges from min to max (options 'off'/'on', 'off' by default)
    67     ParamOut.NbSlice='on'; % edit box nbre of slices made active
     67    ParamOut.NbSlice=1; % impose calculation in a single process (no parallel processing to avoid 'holes'))
    6868    ParamOut.VelType='off';% menu for selecting the velocity type (options 'off'/'one'/'two',  'off' by default)
    6969    ParamOut.FieldName='off';% menu for selecting the field (s) in the input file(options 'off'/'one'/'two', 'off' by default)
    … …  
    8787    if ~exist(FirstFileName,'file')
    8888        msgbox_uvmat('WARNING',['the first input file ' FirstFileName ' does not exist'])
    89     else
    90         [i1,i2,j1,j2] = get_file_index(Param.IndexRange.last_i,last_j,PairString);
    91         LastFileName=fullfile_uvmat(Param.InputTable{1,1},Param.InputTable{1,2},Param.InputTable{1,3},...
    92         Param.InputTable{1,5},Param.InputTable{1,4},i1,i2,j1,j2);
    93         if ~exist(FirstFileName,'file')
    94              msgbox_uvmat('WARNING',['the last input file ' LastFileName ' does not exist'])
    95         end
    9689    end
    9790
    9891    %% check the validity of  input file types
    99     ImageTypeOptions={'image','multimage','mmreader','video'};%allowed input file types(images)
    10092    FileInfo=get_file_info(FirstFileName);
    101     FileType=FileInfo.FileType;
    102     CheckImage=~isempty(find(strcmp(FileType,ImageTypeOptions), 1));% =1 for images
    103     if ~CheckImage
    104         msgbox_uvmat('ERROR',['invalid file type input: ' FileType ' not an image'])
     93    if ~strcmp(FileInfo.FileType,'multimage')
     94        msgbox_uvmat('ERROR',['invalid file type input: ' FileInfo.FileType ' not an image'])
    10595        return
    10696    end
     97    xmlinput=uigetfile_uvmat('pick xml file for timing',fileparts(fileparts(FirstFileName)),'.xml');
     98    [tild,ParamOut.ActionInput.XmlFile]=fileparts(xmlinput);
     99    ParamOut.ActionInput.XmlFile
    107100   
    108101    return
    … …  
    131124ListFile=ListCells(1,find(~check_dir & ~check_bad));
    132125
    133 % %% create the netcdf file with name flname in format NETCDF4
     126%% check file names
     127RootName=regexprep(ListFile{1},'.tif$','')
     128for ilist=2:numel(ListFile)
     129    rank=regexprep(ListFile{ilist},'.tif$','');
     130    rank=regexprep(rank,['^' RootName '@'],'');
     131    if ~isequal(str2num(rank),ilist-1)
     132        disp('error in the list of input files')
     133        return
     134    end
     135end
     136
     137%% output directory
    134138 OutputDir=fullfile(Param.InputTable{1,1},[Param.OutputSubDir Param.OutputDirExt]);
    135139 
    136 % cmode = netcdf.getConstant('NETCDF4');
    137 % cmode = bitor(cmode, netcdf.getConstant('CLASSIC_MODEL'));
    138 % cmode = bitor(cmode, netcdf.getConstant('CLOBBER'));
    139 % nc = netcdf.create(fullfile(OutputDir,'images.nc'), cmode);
     140%% Timing
     141XmlInputFile=fullfile(Param.InputTable{1,1},[Param.ActionInput.XmlFile '.xml'])
     142XmlInput=imadoc2struct(XmlInputFile,'Camera');
    140143
    141144%% Main loop
    142 ImagesPerLevel=100;%100;
     145
     146ImagesPerLevel=size(XmlInput.Time,2)-1;%100;
    143147count=0;
    144148%count=316;%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%CORRECTION EXP08: 4684 images -> start at 316 start 67->_11_1
    … …  
    158162
    159163        if isequal(ImagesPerLevel,1)% mode series
    160            
     164            i_index=count+1;
    161165            OutputFile=fullfile(OutputDir,['img_' num2str(count+1) '.png']);
    162         else % mode multilevel or volume (indices i and j)
     166        else % indices i and j
    163167            i_index=fix(count/ImagesPerLevel)+1;
    164168            j_index=mod(count,ImagesPerLevel)+1;
    … …  
    170174end
    171175
     176%% create the xml file of PCO camera
     177XmlInput.Camera.CameraName='PCO';
     178t=struct2xml(XmlInput.Camera);
     179t=set(t,1,'name','ImaDoc');
     180save(fullfile(Param.InputTable{1,1},'PCO.xml'),t)
    172181
     182%% remove initial files if transfer OK
     183    if i_index== (size(XmlInput.Time,1)-1)
     184
     185        [SUCCESS,MESSAGE]=rmdir(DirImages,'s')
     186       
     187    end

  • trunk/src/series/extract_rdvision.m

    r935 r939  
    338338%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    339339function [BinList,errormsg]=binread_rdv_series(PathDir,SeqData,SqbData,nbfield2,NomTypeNew)
    340 % BINREAD_RDV Permet de lire les fichiers bin gï¿œnï¿œrï¿œs par Hiris ï¿œ partir du
    341 % fichier seq associï¿œ.
     340% BINREAD_RDV Permet de lire les fichiers bin g???n???r???s par Hiris ??? partir du
     341% fichier seq associ???.
    342342%   [IMGS,TIMESTAMPS,NB_FRAMES] = BINREAD_RDV(FILENAME,FRAME_IDX) lit
    343 %   l'image d'indice FRAME_IDX de la sï¿œquence FILENAME.
    344 %
    345 %   Entrï¿œes
     343%   l'image d'indice FRAME_IDX de la s???quence FILENAME.
     344%
     345%   Entr???es
    346346%   -------
    347 %   FILENAME  : Nom du fichier sï¿œquence (.seq).
    348 %   FRAME_IDX : Indice de l'image ï¿œ lire. Si FRAME_IDX vaut -1 alors la
    349 %   sï¿œquence est entiï¿œrement lue. Si FRAME_IDX est un tableau d'indices
     347%   FILENAME  : Nom du fichier s???quence (.seq).
     348%   FRAME_IDX : Indice de l'image ??? lire. Si FRAME_IDX vaut -1 alors la
     349%   s???quence est enti???rement lue. Si FRAME_IDX est un tableau d'indices
    350350%   alors toutes les images d'incides correspondant sont lues. Si FRAME_IDX
    351351%   est un tableau vide alors aucune image n'est lue mais le nombre
    352 %   d'images et tous les timestamps sont renvoyï¿œs. Les indices commencent ï¿œ
    353 %   1 et se termines ï¿œ NB_FRAMES.
     352%   d'images et tous les timestamps sont renvoy???s. Les indices commencent ???
     353%   1 et se termines ??? NB_FRAMES.
    354354%
    355355%   Sorties
    … …  
    357357%   IMGS        : Images de sortie.
    358358%   TIMESTAMPS  : Timestaps des images lues.
    359 %   NB_FRAMES   : Nombres d'images dans la sï¿œquence.
     359%   NB_FRAMES   : Nombres d'images dans la s???quence.
    360360NbBinFile=0;
    361361BinSize=0;
    … …  
    467467end
    468468%% correct NbDti
    469 NbDti=size(timestamp,1); %default for series or burst
     469NbDti=size(timestamp,1)-1; %default for series or burst
    470470uid_motor_nbslice=find(t,'ImaDoc/TranslationMotor/Nbslice');
    471 uid_nbDtk=find(t,'ImaDoc/TranslationMotor/NbDtk');
    472 if ~isempty(uid_motor_nbslice)&& ~isempty(uid_nbDtk)
     471uid_NbDtk=find(t,'ImaDoc/Camera/BurstTiming/NbDtk');
     472if ~isempty(uid_motor_nbslice)&& ~isempty(uid_NbDtk)
    473473    uid_content=get(t,uid_motor_nbslice,'contents');
    474474    NbSlice=str2num(get(t,uid_content,'value'));
    … …  
    496496uid_content_Dtk=get(t,uid_Dtk,'contents');
    497497Dtk=str2num(get(t,uid_content_Dtk,'value'));
    498 uid_NbDtk=find(t,'ImaDoc/Camera/BurstTiming/NbDtk');
    499498uid_content_NbDtk=get(t,uid_NbDtk,'contents');
    500499NbDtk=str2num(get(t,uid_content_NbDtk,'value'));
    … …  
    511510       [success,msg] = fileattrib(newxml,'+w','g');% allow writing access for the group of users 
    512511        if success==0
    513             msgbox_uvmat('WARNING',{['unable to set group write access to ' newxml ':']; msg});%error message for directory creation
     512            disp({['warning: unable to set group write access to ' newxml ':']; msg});%error message for directory creation
    514513        end
    515514

  • trunk/src/series/merge_proj_polar.m

    r938 r939  
    6161%% set the input elements needed on the GUI series when the function is selected in the menu ActionName or InputTable refreshed
    6262if isstruct(Param) && isequal(Param.Action.RUN,0)
    63     ParamOut.AllowInputSort='off';% allow alphabetic sorting of the list of input file SubDir (options 'off'/'on', 'off' by default)
     63    ParamOut.AllowInputSort='on';% allow alphabetic sorting of the list of input file SubDir (options 'off'/'on', 'off' by default)
    6464    ParamOut.WholeIndexRange='off';% prescribes the file index ranges from min to max (options 'off'/'on', 'off' by default)
    6565    ParamOut.NbSlice='off'; %nbre of slices ('off' by default)
    … …  
    6868    ParamOut.FieldTransform = 'on';%can use a transform function
    6969    ParamOut.TransformPath=fullfile(fileparts(which('uvmat')),'transform_field');% path to transform functions (needed for compilation only)
    70     ParamOut.ProjObject='on';%can use projection object(option 'off'/'on',
    71     ParamOut.Mask='on';%can use mask option   (option 'off'/'on', 'off' by default)
    72     ParamOut.OutputDirExt='.mproj';%set the output dir extension
     70    ParamOut.ProjObject='off';%can use projection object(option 'off'/'on',
     71    ParamOut.Mask='off';%can use mask option   (option 'off'/'on', 'off' by default)
     72    ParamOut.OutputDirExt='.polar';%set the output dir extension
    7373    ParamOut.OutputFileMode='NbInput';% '=NbInput': 1 output file per input file index, '=NbInput_i': 1 file per input file index i, '=NbSlice': 1 file per slice
    7474      %check the input files
    … …  
    8383    if ~exist(FirstFileName,'file')
    8484        msgbox_uvmat('WARNING',['the first input file ' FirstFileName ' does not exist'])
    85     elseif isequal(size(Param.InputTable,1),1) && ~isfield(Param,'ProjObject')
    86         msgbox_uvmat('WARNING','You may need a projection object of type plane for merge_proj')
    8785    end
    8886    return
    8987end
     88
     89%%%% specific input parameters
     90% calculate the positions on which to interpolate
     91radius_ref=450;% radius of the mountain top
     92radius_shifted=-120:2:120;% radius shifted by the radius of the origin at the topography summit
     93radius=radius_ref+radius_shifted;%radius from centre of the tank
     94azimuth_arclength=(-120:2:500);%azimuth in arc length at origin position
     95azimuth=pi/2-azimuth_arclength/radius_ref;%azimuth in radian
     96[Radius,Azimuth]=meshgrid(radius,azimuth);
     97XI=Radius.*cos(Azimuth);% set of x axis of the points where interpolqtion needs to be done
     98YI=Radius.*sin(Azimuth)-radius_ref;% set of y axis of the points where interpolqtion needs to be done
     99FieldNames={'vec(U,V)';'curl(U,V)';'div(U,V)'};
     100HeadData.ListVarName= {'radius','azimuth'} ;
     101HeadData.VarDimName={'radius','azimuth'};
     102HeadData.VarAttribute={'coord_y','coord_x'} ;
     103HeadData.radius=radius_shifted;
     104HeadData.azimuth=azimuth_arclength;   
    90105
    91106%%%%%%%%%%%% STANDARD PART (DO NOT EDIT) %%%%%%%%%%%%
    … …  
    105120
    106121%% define the directory for result file (with path=RootPath{1})
    107 OutputDir=['three_cameras' Param.OutputDirExt];% subdirectory for output files
     122OutputDir=[Param.OutputSubDir Param.OutputDirExt];% subdirectory for output files
    108123
    109124if ~isfield(Param,'InputFields')
    … …  
    115130RootFile=Param.InputTable(:,3);
    116131SubDir=Param.InputTable(:,2);
    117 NomType=Param.InputTable(:,4);
     132%NomType=Param.InputTable(:,4);
    118133FileExt=Param.InputTable(:,5);
    119134[filecell,i1_series,i2_series,j1_series,j2_series]=get_file_series(Param);
    … …  
    176191time=mean(time,1); %averaged time taken for the merged field
    177192
     193%% height z
     194    % position of projection plane
     195   
     196ProjObjectCoord=XmlData{1}.GeometryCalib.SliceCoord;
     197CoordUnit=XmlData{1}.GeometryCalib.CoordUnit;
     198for iview =2:numel(XmlData)
     199    if ~(isfield(XmlData{iview},'GeometryCalib')&& isequal(XmlData{iview}.GeometryCalib.SliceCoord,ProjObjectCoord))...
     200        disp('error: geometric calibration missing')
     201        return
     202    end
     203end
     204
     205
    178206%% coordinate transform or other user defined transform
    179207transform_fct='';%default fct handle
    … …  
    194222%% check the validity of  input file types
    195223for iview=1:NbView
    196     if ~isequal(CheckImage{iview},1)&&~isequal(CheckNc{iview},1)
    197         disp_uvmat('ERROR','input set of input series: need  either netcdf either image series',checkrun)
     224    if ~isequal(CheckNc{iview},1)
     225        disp_uvmat('ERROR','input files needs to be in netcdf (extension .nc)',checkrun)
    198226        return
    199227    end
    … …  
    201229
    202230%% output file type
    203 if min(cell2mat(CheckImage))==1 && (~Param.CheckObject || strcmp(Param.ProjObject.Type,'plane'))
    204     FileExtOut='.png'; %image output (input and proj result = image)
    205     for iview=1:NbView
    206         BitDepth(iview)=FileInfo{iview}.BitDepth;
    207     end
    208     BitDepth=max(BitDepth);
    209 else
    210     FileExtOut='.nc'; %netcdf output
    211 end
     231FileExtOut='.nc'; %netcdf output
    212232if isempty(j1_series{1})
    213233    NomTypeOut='_1';
    … …  
    215235    NomTypeOut='_1_1';
    216236end
    217 %NomTypeOut=NomType;% output file index will indicate the first and last ref index in the series
    218237RootFileOut=RootFile{1};
    219238for iview=2:NbView
    … …  
    224243end
    225244
    226 %% mask (TODO: case of multilevels)
    227 MaskData=cell(NbView,1);
    228 if Param.CheckMask
    229     if ischar(Param.MaskTable)% case of a single mask (char chain)
    230         Param.MaskTable={Param.MaskTable};
    231     end
    232     for iview=1:numel(Param.MaskTable)
    233         if exist(Param.MaskTable{iview},'file')
    234             [MaskData{iview},tild,errormsg] = read_field(Param.MaskTable{iview},'image');
    235             if ~isempty(transform_fct) && nargin(transform_fct)>=2
    236                 MaskData{iview}=transform_fct(MaskData{iview},XmlData{iview});
    237             end
    238         end
    239     end
    240 end
    241 
    242 %% Set field names and velocity types
    243 %use Param.InputFields for all views
    244245
    245246%% MAIN LOOP ON FIELDS
    … …  
    250251%     nbmissing=0;
    251252
    252     %%%%%%%%%%%%%%%% loop on field indices %%%%%%%%%%%%%%%%
     253%%%%%%%%%%%%%%%% loop on field indices %%%%%%%%%%%%%%%%
    253254tstart=tic; %used to record the computing time
    254255CheckOverwrite=1;%default
    … …  
    257258end
    258259for index=1:NbField
    259         update_waitbar(WaitbarHandle,index/NbField)
     260    update_waitbar(WaitbarHandle,index/NbField)
    260261    if ~isempty(RUNHandle) && ~strcmp(get(RUNHandle,'BusyAction'),'queue')
    261262        disp('program stopped by user')
    … …  
    263264    end
    264265   
    265      %% generating the name of the merged field
     266    %% generating the name of the merged field
    266267    i1=i1_series{1}(index);
    267268    if ~isempty(i2_series{end})
    … …  
    282283    OutputFile=fullfile_uvmat(RootPath{1},OutputDir,RootFileOut,FileExtOut,NomTypeOut,i1,i2,j1,j2);
    283284    if ~CheckOverwrite && exist(OutputFile,'file')
    284             disp(['existing output file ' OutputFile ' already exists, skip to next field'])
    285             continue% skip iteration if the mode overwrite is desactivated and the result file already exists
    286     end
    287        
     285        disp(['existing output file ' OutputFile ' already exists, skip to next field'])
     286        continue% skip iteration if the mode overwrite is desactivated and the result file already exists
     287    end
     288   
     289    %% z position
     290       ZIndex=mod(i1_series{1}(index)-1,NbSlice_calib{1})+1;%Zindex for phys transform
     291       ZPos=ProjObjectCoord(ZIndex,3);
     292       % radius of the topography section at z pos
     293       TopoRadius=0;
     294       ind_mask=[];
     295       if ZPos<20
     296    TopoRadius=40*sin(acos((20+ZPos)/40));
     297       ind_mask=(XI'.*XI'+YI'.*YI')<TopoRadius*TopoRadius;% indidces of data to mask
     298       end
     299       
    288300    %%%%%%%%%%%%%%%% loop on views (input lines) %%%%%%%%%%%%%%%%
    289301    Data=cell(1,NbView);%initiate the set Data
    290302    timeread=zeros(1,NbView);
    291303    for iview=1:NbView
    292         %% reading input file(s)     
     304        %% reading input file(s)
    293305        [Data{iview},tild,errormsg] = read_field(filecell{iview,index},FileType{iview},ParamIn{iview},frame_index{iview}(index));
    294306        if ~isempty(errormsg)
    … …  
    318330        end
    319331       
    320         %% calculate tps coefficients if needed
    321         check_proj_tps= 1;
    322         Data{iview}=tps_coeff_field(Data{iview},check_proj_tps);
     332        %% calculate tps coefficients
     333        Data{iview}=tps_coeff_field(Data{iview},1);
    323334       
    324         %% projection on object (gridded plane)
    325         if Param.CheckObject
    326             [Data{iview},errormsg]=proj_field(Data{iview},Param.ProjObject);
    327             if ~isempty(errormsg)
    328                 disp_uvmat('ERROR',['ERROR in merge_proge/proj_field: ' errormsg],checkrun)
    329                 return
    330             end
    331         end
    332        
     335        %% projection on the polar grid
     336        [DataOut,VarAttribute,errormsg]=calc_field_tps(Data{iview}.Coord_tps,Data{iview}.NbCentre,Data{iview}.SubRange,...
     337            cat(3,Data{iview}.U_tps,Data{iview}.V_tps),FieldNames,cat(3,XI,YI));
     338        ListVarName=(fieldnames(DataOut))';
     339        ProjData{iview}=HeadData;
     340        ProjData{iview}.ListVarName= [ProjData{iview}.ListVarName ListVarName];
     341        ProjData{iview}.VarDimName={'radius','azimuth'};
     342        ProjData{iview}.VarAttribute=[{'coord_x'} {'coord_y'} VarAttribute];
     343        for ivar=1:numel(ListVarName)
     344            ProjData{iview}.VarDimName{ivar+2}={'radius','azimuth'};
     345            ProjData{iview}.(ListVarName{ivar})=(DataOut.(ListVarName{ivar}))';
     346        end
    333347        %% mask
    334         if Param.CheckMask && ~isempty(MaskData{iview})
    335              [Data{iview},errormsg]=mask_proj(Data{iview},MaskData{iview});
    336         end
     348        %         if Param.CheckMask && ~isempty(MaskData{iview})
     349        %              [Data{iview},errormsg]=mask_proj(Data{iview},MaskData{iview});
     350        %         end
    337351    end
    338352    %%%%%%%%%%%%%%%% END LOOP ON VIEWS %%%%%%%%%%%%%%%%
    339 
     353   
    340354    %% merge the NbView fields
    341     [MergeData,errormsg]=merge_field(Data);
     355    [MergeData,errormsg]=merge_field(ProjData);
    342356    if ~isempty(errormsg)
    343357        disp_uvmat('ERROR',errormsg,checkrun);
    344358        return
    345359    end
    346 
     360   
    347361    %% time of the merged field: take the average of the different views
    348362    if ~isempty(time)
    349         timeread=time(index);   
     363        timeread=time(index);
    350364    elseif ~isempty(find(timeread))% time defined from ImaDoc
    351365        timeread=mean(timeread(timeread~=0));% take average over times form the files (when defined)
    352366    else
    353         timeread=index;% take time=file index
    354     end
    355 
    356  
     367        timeread=index;% take time=file index
     368    end
     369   
     370    %% rotating the velocity vectors to the local axis of the polatr coordinates
     371    Unew=MergeData.U.*sin(Azimuth')-MergeData.V.*cos(Azimuth');
     372    MergeData.V=MergeData.U.*cos(Azimuth')+MergeData.V.*sin(Azimuth');
     373    MergeData.U=Unew;
     374    if ~isempty(ind_mask)
     375        MergeData.U(ind_mask)=NaN;
     376        MergeData.V(ind_mask)=NaN;
     377        MergeData.curl(ind_mask)=NaN;
     378        MergeData.div(ind_mask)=NaN;
     379    end
     380   
    357381    %% recording the merged field
    358     if strcmp(FileExtOut,'.png')    %output as image
    359         if BitDepth==8
    360             imwrite(uint8(MergeData.A),OutputFile,'BitDepth',8)
    361         else
    362             imwrite(uint16(MergeData.A),OutputFile,'BitDepth',16)
    363         end
    364         if index==1
    365             %write xml calibration file, using the first file
    366             siz=size(MergeData.A);
    367             npy=siz(1);
    368             npx=siz(2);
    369             if isfield(MergeData,'Coord_x') && isfield(MergeData,'Coord_y')
    370                 Rangx=MergeData.Coord_x;
    371                 Rangy=MergeData.Coord_y;
    372             elseif isfield(MergeData,'AX')&& isfield(MergeData,'AY')
    373                 Rangx=[MergeData.AX(1) MergeData.AX(end)];
    374                 Rangy=[MergeData.AY(1) MergeData.AY(end)];
    375             else
    376                 Rangx=[0.5 npx-0.5];
    377                 Rangy=[npy-0.5 0.5];%default
     382   
     383    MergeData.ListGlobalAttribute={'Conventions','Project','InputFile_1','InputFile_end','NbCoord','NbDim','CoordUnit','ZPos'};
     384    MergeData.Conventions='uvmat';
     385    MergeData.NbCoord=2;
     386    MergeData.NbDim=2;
     387    MergeData.CoordUnit=CoordUnit;
     388    MergeData.ZPos=ZPos;
     389    % time interval of PIV
     390    Dt=[];
     391    if isfield(Data{1},'Dt')&& isnumeric(Data{1}.Dt)
     392        Dt=Data{1}.Dt;
     393    end
     394    for iview =2:numel(Data)
     395        if ~(isfield(Data{iview},'Dt')&& isequal(Data{iview}.Dt,Dt))
     396            Dt=[];%dt not the same for all fields
     397        end
     398    end
     399    if ~isempty(timeread)
     400        MergeData.ListGlobalAttribute=[MergeData.ListGlobalAttribute {'Time'}];
     401        MergeData.Time=timeread;
     402    end
     403
     404    % time unit
     405    if isfield(Data{1},'TimeUnit')
     406        TimeUnit=Data{1}.TimeUnit;
     407        for iview =2:numel(Data)
     408            if ~(isfield(Data{iview},'TimeUnit')&& isequal(Data{iview}.TimeUnit,TimeUnit))
     409                TimeUnit=[];%TimeUnit not the same for all fields
    378410            end
    379             pxcmx=(npx-1)/(Rangx(2)-Rangx(1));
    380             pxcmy=(npy-1)/(Rangy(1)-Rangy(2));
    381             T_x=-pxcmx*Rangx(1)+0.5;
    382             T_y=-pxcmy*Rangy(2)+0.5;
    383             GeometryCal.CalibrationType='rescale';
    384             GeometryCal.CoordUnit=MergeData.CoordUnit;
    385             GeometryCal.focal=1;
    386             GeometryCal.R=[pxcmx,0,0;0,pxcmy,0;0,0,1];
    387             GeometryCal.Tx_Ty_Tz=[T_x T_y 1];
    388             ImaDoc.GeometryCalib=GeometryCal;
    389             t=struct2xml(ImaDoc);
    390             t=set(t,1,'name','ImaDoc');
    391             save(t,[fileparts(OutputFile) '.xml'])
    392         end
    393        
    394     else   %output as netcdf files
    395         MergeData.ListGlobalAttribute={'Conventions','Project','InputFile_1','InputFile_end','NbCoord','NbDim'};
    396         MergeData.Conventions='uvmat';
    397         MergeData.NbCoord=2;
    398         MergeData.NbDim=2;
    399         % time interval of PIV
    400         Dt=[];
    401         if isfield(Data{1},'Dt')&& isnumeric(Data{1}.Dt)
    402             Dt=Data{1}.Dt;
    403         end
    404         for iview =2:numel(Data)
    405             if ~(isfield(Data{iview},'Dt')&& isequal(Data{iview}.Dt,Dt))
    406                 Dt=[];%dt not the same for all fields
    407             end
    408         end
    409         if ~isempty(timeread)
    410             MergeData.ListGlobalAttribute=[MergeData.ListGlobalAttribute {'Time'}];
    411             MergeData.Time=timeread;
    412         end
    413         % position of projection plane
    414         if isfield(Data{1},'ProjObjectCoord')&& isfield(Data{1},'ProjObjectAngle')
    415             'test'
    416             ProjObjectCoord=Data{1}.ProjObjectCoord;
    417             ProjObjectAngle=Data{1}.ProjObjectAngle;
    418             for iview =2:numel(Data)
    419                 if ~(isfield(Data{iview},'ProjObjectCoord')&& isequal(Data{iview}.ProjObjectCoord,ProjObjectCoord))...
    420                         ||~(isfield(Data{iview},'ProjObjectAngle')&& isequal(Data{iview}.ProjObjectAngle,ProjObjectAngle))
    421                     ProjObjectCoord=[];%dt not the same for all fields
    422                 end
    423             end
    424             if ~isempty(ProjObjectCoord)
    425                 MergeData.ListGlobalAttribute=[MergeData.ListGlobalAttribute {'ProjObjectCoord'} {'ProjObjectAngle'}];
    426                 MergeData.ProjObjectCoord=ProjObjectCoord;
    427                 MergeData.ProjObjectAngle=ProjObjectAngle;
    428             end
    429         end
    430         % coord unit
    431         if isfield(Data{1},'CoordUnit')
    432             CoordUnit=Data{1}.CoordUnit;
    433             for iview =2:numel(Data)
    434                 if ~(isfield(Data{iview},'CoordUnit')&& isequal(Data{iview}.CoordUnit,CoordUnit))
    435                     CoordUnit=[];%CoordUnit not the same for all fields
    436                 end
    437             end
    438             if ~isempty(CoordUnit)
    439                 MergeData.ListGlobalAttribute=[MergeData.ListGlobalAttribute {'CoordUnit'}];
    440                 MergeData.CoordUnit=CoordUnit;
    441             end
    442         end
    443         % time unit
    444         if isfield(Data{1},'TimeUnit')
    445             TimeUnit=Data{1}.TimeUnit;
    446             for iview =2:numel(Data)
    447                 if ~(isfield(Data{iview},'TimeUnit')&& isequal(Data{iview}.TimeUnit,TimeUnit))
    448                     TimeUnit=[];%TimeUnit not the same for all fields
    449                 end
    450             end
    451             if ~isempty(TimeUnit)
    452                 MergeData.ListGlobalAttribute=[MergeData.ListGlobalAttribute {'TimeUnit'}];
    453                 MergeData.TimeUnit=TimeUnit;
    454             end
    455         end
    456         error=struct2nc(OutputFile,MergeData);%save result file
    457         if isempty(error)
    458             disp(['output file ' OutputFile ' written'])
    459         else
    460             disp(error)
    461         end
    462     end
    463 end
     411        end
     412        if ~isempty(TimeUnit)
     413            MergeData.ListGlobalAttribute=[MergeData.ListGlobalAttribute {'TimeUnit'}];
     414            MergeData.TimeUnit=TimeUnit;
     415        end
     416    end
     417    error=struct2nc(OutputFile,MergeData);%save result file
     418    if isempty(error)
     419        disp(['output file ' OutputFile ' written'])
     420    else
     421        disp(error)
     422    end
     423end
     424
    464425ellapsed_time=toc(tstart);
    465426disp(['total ellapsed time ' num2str(ellapsed_time/60,2) ' minutes'])

  • trunk/src/struct2xml.m

    r924 r939  
    1111% OUTPUT:
    1212% t: xmltree reproducing the structure of Object
    13 % type 'save(t)' to visualize the xml text and save(filename,t) to save it in a file
     13% type 'save(t)' to visualize the xml text and save(t,filename) to save it in a file
    1414% to set the title tag of the file, type t=set(t,1,'name','titletag')
    1515%
Note: See TracChangeset for help on using the changeset viewer.