Changeset 905 for trunk/src


Ignore:
Timestamp:
May 30, 2015, 8:44:20 PM (9 years ago)
Author:
sommeria
Message:

projection with tps corrected + minor changes

Location:
trunk/src
Files:
7 edited

Legend:

Unmodified
Added
Removed
  • trunk/src/calc_field_tps.m

    r867 r905  
    113113    nbvec_sub=NbCentre(isub);
    114114    check_range=(Coord_interp >=ones(nb_sites,1)*SubRange(:,1,isub)' & Coord_interp<=ones(nb_sites,1)*SubRange(:,2,isub)');
    115     ind_sel=find(sum(check_range,2)==nb_coord);
     115    ind_sel=find(sum(check_range,2)==nb_coord);% select points whose all coordinates are in the prescribed range
    116116    nbval(ind_sel)=nbval(ind_sel)+1;% records the number of values for eacn interpolation point (in case of subdomain overlap)
    117117    if check_grid
  • trunk/src/mouse_up.m

    r871 r905  
    141141        if  ~isempty(ObjectData)
    142142            % plot the field projected on the object
    143             ProjData= proj_field(UvData.Field,ObjectData);%project the current interface field on ObjectData
    144             if ~isempty(ProjData)
     143            [ProjData,errormsg]= proj_field(UvData.Field,ObjectData);%project the current interface field on ObjectData
     144            if isempty(errormsg) && ~isempty(ProjData)
    145145                if strcmp(FigTag,'uvmat')% uvmat plot selected, projection plot seen in view_field
    146146                    hview_field=findobj(allchild(0),'tag','view_field');
     
    170170                end
    171171            end
     172            if ~isempty(errormsg)
     173                msgbox_uvmat('ERROR',errormsg)
     174                return
     175            end
    172176            set(hhuvmat.CheckViewField,'Value',1);%
    173177            set(hhuvmat.CheckEditObject,'Value',1);%   
  • trunk/src/proj_field.m

    r896 r905  
    522522function  [ProjData,errormsg] = proj_line(FieldData, ObjectData)
    523523%-----------------------------------------------------------------
     524
     525%% prepare heading for the projected field
    524526[ProjData,errormsg]=proj_heading(FieldData,ObjectData);%transfer global attributes
    525527if ~isempty(errormsg)
     
    544546ListIndex={};
    545547
     548%% group the variables (fields of 'FieldData') in cells of variables with the same dimensions
     549[CellInfo,NbDim,errormsg]=find_field_cells(FieldData);
     550if ~isempty(errormsg)
     551    errormsg=['error in proj_field/proj_line:' errormsg];
     552    return
     553end
     554CellInfo=CellInfo(NbDim==2); %keep only the 2D cells
     555%%%%%% TODO: treat 1D fields: project as identity so that P o P=P for projection operation
     556cell_select=true(size(CellInfo));
     557
     558for icell=1:length(CellInfo)
     559    if isfield(CellInfo{icell},'ProjModeRequest')
     560        if ~strcmp(CellInfo{icell}.ProjModeRequest, ProjMode)
     561            cell_select(icell)=0;
     562        end
     563        if strcmp(ProjMode,'interp_tps')&& ~strcmp(CellInfo{icell}.CoordType,'tps')
     564            cell_select(icell)=0;
     565        end
     566    end
     567end
     568if isempty(find(cell_select))
     569    errormsg=[' invalid projection mode ''' ProjMode ''': use ''interp_tps'' to interpolate spatial derivatives'];
     570    return
     571end
     572CellInfo=CellInfo(cell_select);
    546573
    547574%% projection line: object types selected from  proj_field='line','polyline','polygon','rectangle','ellipse':
     
    633660    end
    634661end
    635 
    636 
    637 %% group the variables (fields of 'FieldData') in cells of variables with the same dimensions
    638 [CellInfo,NbDim,errormsg]=find_field_cells(FieldData);
    639 if ~isempty(errormsg)
    640     errormsg=['error in proj_field/proj_line:' errormsg];
    641     return
    642 end
    643 CellInfo=CellInfo(NbDim==2); %keep only the 2D cells
    644 %%%%%% TODO: treat 1D fields: project as identity so that P o P=P for projection operation
    645662
    646663%% loop on variable cells with the same space dimension 2
     
    12871304                    Dist=Distx.*Distx+Disty.*Disty;
    12881305                    for ivar=1:numel(VarVal)
    1289                         VarVal{ivar}(Dist>4*ProjData.CoordMesh)=NaN;% put to NaN interpolated positions too far from initial data
     1306                        VarVal{ivar}(Dist>16*ProjData.CoordMesh)=NaN;% % put to NaN interpolated positions further than 4 meshes from initial data
    12901307                    end 
    12911308                   
     
    13251342               
    13261343                % set to NaN interpolation points which are too far from any initial data (more than 2 CoordMesh)
    1327                     if exist('scatteredInterpolant','file')%recent Matlab versions
    1328                         F=scatteredInterpolant(coord_X, coord_Y,coord_X,'nearest');
    1329                         G=scatteredInterpolant(coord_X, coord_Y,coord_Y,'nearest');
    1330                     else
    1331                         F=TriScatteredInterp([coord_X coord_Y],coord_X,'nearest');
    1332                         G=TriScatteredInterp([coord_X coord_Y],coord_Y,'nearest');
    1333                     end
    1334                     Distx=F(XI,YI)-XI;% diff of x coordinates with the nearest measurement point
    1335                     Disty=G(XI,YI)-YI;% diff of y coordinates with the nearest measurement point
    1336                     Dist=Distx.*Distx+Disty.*Disty;
    1337                     for ivar=1:numel(VarVal)
    1338                         VarVal{ivar}(Dist>4*ProjData.CoordMesh)=NaN;% put to NaN interpolated positions too far from initial data
    1339                     end 
    1340                    
    1341                
     1344                Coord=permute(Coord,[1 3 2]);
     1345                Coord=reshape(Coord,size(Coord,1)*size(Coord,2),2);
     1346                if exist('scatteredInterpolant','file')%recent Matlab versions
     1347                    F=scatteredInterpolant(Coord,Coord(:,1),'nearest');
     1348                    G=scatteredInterpolant(Coord,Coord(:,2),'nearest');
     1349                else
     1350                    F=TriScatteredInterp(Coord,Coord(:,1),'nearest');
     1351                    G=TriScatteredInterp(Coord,Coord(:,2),'nearest');
     1352                end
     1353                Distx=F(XI,YI)-XI;% diff of x coordinates with the nearest measurement point
     1354                Disty=G(XI,YI)-YI;% diff of y coordinates with the nearest measurement point
     1355                Dist=Distx.*Distx+Disty.*Disty;
    13421356                ListVarName=(fieldnames(DataOut))';
    13431357                VarDimName=cell(size(ListVarName));
     
    13451359                    VarName=ListVarName{ilist};
    13461360                    ProjData.(VarName)=DataOut.(VarName);
     1361                    ProjData.(VarName)(Dist>16*ProjData.CoordMesh)=NaN;% put to NaN interpolated positions further than 4 meshes from initial data
    13471362                    VarDimName{ilist}={'coord_y','coord_x'};
    13481363                end
  • trunk/src/series.m

    r904 r905  
    16071607    detect=exist(fullfile(Param.InputTable{1,1},SubDirOutNew),'dir');% test if  the dir  already exist
    16081608    check_create=1; %need to create the result directory by default
     1609    CheckOverwrite=1;
     1610    if isfield(Param,'CheckOverwrite')
     1611        CheckOverwrite=Param.CheckOverwrite;
     1612    end
    16091613    while detect
    1610         if Param.CheckOverwrite
     1614        if CheckOverwrite
    16111615            comment=', possibly overwrite previous data';
    16121616        else
     
    17191723end
    17201724nbfield_j=numel(ref_j); % number of j indices
     1725BlockLength=numel(ref_i);%default
    17211726if isempty(Param.IndexRange.NbSlice)
    17221727    NbProcess=NbCore;% choose one process per core by default if NbSlice is not imposed
     
    18051810    %% processing on a different session of the same computer (background) or cluster, create executable files
    18061811    batch_file_list=cell(NbProcess,1);% initiate the list of executable files
    1807     DirBat=fullfile(OutputDir,'0_EXE');
     1812    DirExe=fullfile(OutputDir,'0_EXE');%directory name for executable files
    18081813    switch computer
    18091814        case {'PCWIN','PCWIN64'} %Windows system
     
    18131818    end
    18141819    %create subdirectory for executable files
    1815     if ~exist(DirBat,'dir')
    1816         [tild,msg1]=mkdir(DirBat);
     1820    if ~exist(DirExe,'dir')
     1821        [tild,msg1]=mkdir(DirExe);
    18171822        if ~strcmp(msg1,'')
    1818             errormsg=['cannot create ' DirBat ': ' msg1];%error message for directory creation
     1823            errormsg=['cannot create ' DirExe ': ' msg1];%error message for directory creation
    18191824            return
    18201825        end
     
    18521857       
    18531858        %create the executable file
    1854          filebat=fullfile_uvmat(DirBat,'',Param.InputTable{1,3},ExeExt,OutputNomType,...
     1859         filebat=fullfile_uvmat(DirExe,'',Param.InputTable{1,3},ExeExt,OutputNomType,...
    18551860           Param.IndexRange.first_i,Param.IndexRange.last_i,first_j,last_j);
    18561861        batch_file_list{iprocess}=filebat;
     
    19191924    case 'cluster_oar' % option 'oar-parexec' used
    19201925        %create subdirectory for oar command and log files
    1921         DirOAR=fullfile(OutputDir,'0_LOG');
    1922         if exist(DirOAR,'dir')% delete the content of the dir 0_LOG to allow new input
    1923             curdir=pwd;
    1924             cd(DirOAR)
    1925             delete('*')
    1926             cd(curdir)
    1927         else
    1928             [tild,msg1]=mkdir(DirOAR);
    1929             if ~strcmp(msg1,'')
    1930                 errormsg=['cannot create ' DirOAR ': ' msg1];%error message for directory creation
    1931                 return
    1932             end
    1933         end
    1934         filename_joblist=fullfile(DirOAR,'0_job_list.txt');%create name of the global executable file
     1926        %DirOARLog=fullfile(OutputDir,'0_LOG');
     1927        %DirOARExe=fullfile(OutputDir,'0_EXE');
     1928%         if exist(DirOAR,'dir')% delete the content of the dir 0_LOG to allow new input
     1929%             curdir=pwd;
     1930%             cd(DirOAR)
     1931%             delete('*')
     1932%             cd(curdir)
     1933%         else
     1934%             [tild,msg1]=mkdir(DirOAR);
     1935%             if ~strcmp(msg1,'')
     1936%                 errormsg=['cannot create ' DirOAR ': ' msg1];%error message for directory creation
     1937%                 return
     1938%             end
     1939%         end
     1940        filename_joblist=fullfile(DirExe,'0_job_list.txt');%create name of the global executable file
     1941        filename_log=fullfile(DirLog,'0_job_list.stdout');%file for output messages of the master oar process
     1942        filename_errors=fullfile(DirLog,'0_job_list.stderr');%file for error messages of the master oar process
     1943       
    19351944        fid=fopen(filename_joblist,'w');
    19361945        for p=1:length(batch_file_list)
     
    19581967            '-l /core=' num2str(NbCore) ','...
    19591968            'walltime=' datestr(WallTimeOneJob/24,13) ' '...
    1960             '-E ' regexprep(filename_joblist,'\.txt\>','.stderr') ' '...
    1961             '-O ' regexprep(filename_joblist,'\.txt\>','.stdout') ' '...
     1969            '-E ' filename_errors ' '...
     1970            '-O ' filename_log ' '...
    19621971            extra_oar ' '...
    19631972            '"oar-parexec -s -f ' filename_joblist ' '...
    19641973            '-l ' filename_joblist '.log"\n'];
    19651974       
    1966         filename_oarcommand=fullfile(DirOAR,'0_oar_command');
     1975        filename_oarcommand=fullfile(DirExe,'0_oar_command');
    19671976        fid=fopen(filename_oarcommand,'w');
    19681977        fprintf(fid,oar_command);
     
    35973606% --- Executes on button press in CheckOverwrite.
    35983607function CheckOverwrite_Callback(hObject, eventdata, handles)
    3599 % hObject    handle to CheckOverwrite (see GCBO)
    3600 % eventdata  reserved - to be defined in a future version of MATLAB
    3601 % handles    structure with handles and user data (see GUIDATA)
    3602 
    3603 % Hint: get(hObject,'Value') returns toggle state of CheckOverwrite
     3608
    36043609
    36053610
  • trunk/src/series/civ_series.m

    r904 r905  
    333333    if ~isempty(RUNHandle)% update the waitbar in interactive mode with GUI series  (checkrun=1)
    334334        update_waitbar(WaitbarHandle,ifield/NbField)
    335         if  ~strcmp(get(RUNHandle,'BusyAction'),'queue')
     335        if  checkrun && ~strcmp(get(RUNHandle,'BusyAction'),'queue')
    336336            disp('program stopped by user')
    337337            break
  • trunk/src/series/merge_proj.m

    r883 r905  
    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
     75    ParamOut.CheckOverwriteVisible='on'; % manage the overwrite of existing files (default=1)
    7576    first_j=[];
    7677    if isfield(Param.IndexRange,'first_j'); first_j=Param.IndexRange.first_j; end
     
    9091%%%%%%%%%%%% STANDARD PART (DO NOT EDIT) %%%%%%%%%%%%
    9192ParamOut=[]; %default output
     93RUNHandle=[];
     94WaitbarHandle=[];
    9295%% read input parameters from an xml file if input is a file name (batch mode)
    9396checkrun=1;
     
    9598    Param=xml2struct(Param);% read Param as input file (batch case)
    9699    checkrun=0;
    97 end
    98 hseries=findobj(allchild(0),'Tag','series');
    99 RUNHandle=findobj(hseries,'Tag','RUN');%handle of RUN button in GUI series
    100 WaitbarHandle=findobj(hseries,'Tag','Waitbar');%handle of waitbar in GUI series
     100else
     101    hseries=findobj(allchild(0),'Tag','series');
     102    RUNHandle=findobj(hseries,'Tag','RUN');%handle of RUN button in GUI series
     103    WaitbarHandle=findobj(hseries,'Tag','Waitbar');%handle of waitbar in GUI series
     104end
    101105
    102106%% define the directory for result file (with path=RootPath{1})
     
    407411    end
    408412end
    409 disp(['total ellapsed time ' num2str(toc(tstart))])
     413ellapsed_time=toc(tstart);
     414disp(['total ellapsed time ' num2str(ellapsed_time/60,2) ' minutes'])
     415disp([ num2str(ellapsed_time/(60*NbField),3) ' minutes per iteration'])
    410416
    411417%'merge_field': concatene fields
  • trunk/src/uvmat.m

    r894 r905  
    37353735        [ObjectData,errormsg]=proj_field(UvData.Field,UvData.ProjObject{iobj});% project field on the object
    37363736        if ~isempty(errormsg)
     3737            errormsg=['projection on ' UvData.ProjObject{iobj}.Type ': ' errormsg ];
    37373738            return
    37383739        end
     
    53785379
    53795380%% update the projection plot on uvmat
    5380 ProjData= proj_field(UvData.Field,ObjectData);%project the current input field on object ObjectData
     5381[ProjData,errormsg]= proj_field(UvData.Field,ObjectData);%project the current input field on object ObjectData
     5382    if ~isempty(errormsg)
     5383        msgbox_uvmat('ERROR',['projection on ' UvData.ProjObject{iobj}.Type ': ' errormsg ])
     5384        return
     5385    end
    53815386plot_field(ProjData,handles.PlotAxes,read_GUI(handles.uvmat));% plot the projected field;
    53825387UvData.PlotAxes=ProjData;% store the plotted field for further update
     
    54465451    Data=get(hview_field,'UserData');
    54475452    hhview_field=guidata(hview_field);
    5448     ProjData= proj_field(UvData.Field,ObjectData);%project the current interface field on ObjectData
     5453    [ProjData,errormsg]= proj_field(UvData.Field,ObjectData);%project the current interface field on ObjectData
     5454    if ~isempty(errormsg)
     5455        msgbox_uvmat('ERROR',['projection on ' UvData.ProjObject{iobj}.Type ': ' errormsg ])
     5456        return
     5457    end
    54495458    [PlotType,PlotParam]=plot_field(ProjData,hhview_field.PlotAxes,read_GUI(hview_field));%read plotting parameters on the uvmat interface
    54505459    haxes=findobj(hview_field,'tag','axes3');
     
    56145623   
    56155624    %% show the projection of the selected object on view_field
    5616     ProjData= proj_field(UvData.Field,UvData.ProjObject{IndexObj});%project the current field on ObjectData
     5625    [ProjData,errormsg]= proj_field(UvData.Field,UvData.ProjObject{IndexObj});%project the current field on ObjectData
     5626    if ~isempty(errormsg)
     5627        msgbox_uvmat('ERROR',['projection on ' UvData.ProjObject{IndexObj}.Type ': ' errormsg])
     5628        return
     5629    end
    56175630    hview_field=findobj(allchild(0),'tag','view_field');
    56185631    if isempty(hview_field)
Note: See TracChangeset for help on using the changeset viewer.