[646]  1  %'mask_proj': restrict input fields to a mask region, set to 0 outside


[639]  2  %


[646]  3  % function [ProjData,errormsg]=mask_proj(FieldData,MaskData)


[639]  4  %


 5  %AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA


 6  % Copyright Joel Sommeria, 2008, LEGI / CNRSUJFINPG, sommeria@coriolislegi.org.


 7  %AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA


 8  % This file is part of the toolbox UVMAT.


 9  %


 10  % UVMAT is free software; you can redistribute it and/or modify


 11  % it under the terms of the GNU General Public License as published by


 12  % the Free Software Foundation; either version 2 of the License, or


 13  % (at your option) any later version.


 14  %


 15  % UVMAT is distributed in the hope that it will be useful,


 16  % but WITHOUT ANY WARRANTY; without even the implied warranty of


 17  % MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the


 18  % GNU General Public License (file UVMAT/COPYING.txt) for more details.


 19  %AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAA


 20 


 21  function [ProjData,errormsg]=mask_proj(FieldData,MaskData)


 22  errormsg='';%default


 23  ProjData=FieldData;


 24 


 25 


 26  %% group the variables (fields of 'FieldData') in cells of variables with the same dimensions


 27  [CellInfo,NbDimArray,errormsg]=find_field_cells(FieldData);


 28  if ~isempty(errormsg)


 29  errormsg=['error in proj_field/proj_plane:' errormsg];


 30  return


 31  end


 32  [Npy,Npx]=size(MaskData.A);


[675]  33 


[639]  34  for icell=1:numel(CellInfo)


 35  if NbDimArray(icell)==2


[672]  36  if isfield(CellInfo{icell},'VarIndex_errorflag')


 37  FFName=FieldData.ListVarName{CellInfo{icell}.VarIndex_errorflag};


[675]  38  check_new_error_flag=0;


[672]  39  else


 40  FFName='FF';%default error (mask) flag name (if not already used)


 41  if isfield(FieldData,'FF')


 42  ind=1;


 43  while isfield(FieldData,['FF_' num2str(ind)])


 44  ind=ind+1;


 45  end


 46  FFName=['FF_' num2str(ind)];% append an index to the name of error flag, FF_1,FF_2...


 47  end


 48  ProjData.ListVarName=[FieldData.ListVarName {FFName}];


[675]  49  ProjData.VarAttribute{numel(ProjData.ListVarName)}.Role='errorflag';


 50  check_new_error_flag=1;


[672]  51  end


[639]  52  switch CellInfo{icell}.CoordType;


[675]  53  case 'scattered'


 54  XName=FieldData.ListVarName{CellInfo{icell}.Coord_x};


 55  YName=FieldData.ListVarName{CellInfo{icell}.Coord_y};


[782]  56  DX=(MaskData.Coord_x(2)MaskData.Coord_x(1))/(Npx1);


 57  DY=(MaskData.Coord_y(2)MaskData.Coord_y(1))/(Npy1);


 58  mask_ind_i=round(0.5+(FieldData.(XName)MaskData.Coord_x(1))/DX);%nbpoint terms


 59  mask_ind_j=round(0.5+(FieldData.(YName)MaskData.Coord_y(1))/DY);%nbpoint terms


[639]  60  checkin=mask_ind_j+Npy*(mask_ind_i1);%array of mask indices for the nbpoints


 61  checkin=checkin(mask_ind_i>=1 & mask_ind_i<=Npx & mask_ind_j>=1 & mask_ind_j<=Npy);%reduced array of mask indices (inside the image)


 62  checkfalse=true(size(FieldData.(XName)));


 63  MaskData.A=reshape(MaskData.A,1,[]);


 64  checkfalse(MaskData.A(checkin)>200)=0;


 65  for ivar=1:numel(CellInfo{icell}.VarIndex)


 66  VarName=FieldData.ListVarName{CellInfo{icell}.VarIndex(ivar)};


[672]  67  ProjData.(VarName)(checkfalse)=0;


[639]  68  end


[675]  69  if check_new_error_flag% an error flag needs to be created in the current cell


[672]  70  ProjData.(FFName)=zeros(size(ProjData.(VarName)));


[675]  71  ProjData.VarDimName=[FieldData.VarDimName FieldData.VarDimName(CellInfo{icell}.VarIndex(1))];


[672]  72  end


[675]  73  ProjData.(FFName)(checkfalse)=1;% update the existing error flag


[639]  74  case 'grid'


[675]  75  XName=FieldData.ListVarName{CellInfo{icell}.CoordIndex(2)};


 76  YName=FieldData.ListVarName{CellInfo{icell}.CoordIndex(1)};


[639]  77  Var1Name=FieldData.ListVarName{CellInfo{icell}.VarIndex(1)};


 78  [Npy_field,Npx_field]=size(FieldData.(Var1Name));


[675]  79  XArray=linspace(FieldData.(XName)(1),FieldData.(XName)(end),Npx_field);


 80  YArray=linspace(FieldData.(YName)(1),FieldData.(YName)(end),Npy_field);


[782]  81  XMask=linspace(MaskData.Coord_x(1),MaskData.Coord_x(end),Npx);


 82  YMask=linspace(MaskData.Coord_y(1),MaskData.Coord_y(end),Npy);


[639]  83  [XMask,YMask]=meshgrid(XMask,YMask);


 84  Mask = interp2(XMask,YMask,MaskData.A,XArray,YArray','nearest');


[672]  85  Mask=Mask>200;


[639]  86  for ivar=1:numel(CellInfo{icell}.VarIndex)


 87  VarName=FieldData.ListVarName{CellInfo{icell}.VarIndex(ivar)};


[672]  88  if ~strcmp(VarName,FFName)


 89  ProjData.(VarName)=FieldData.(VarName).*Mask;


 90  end


[639]  91  end


[675]  92  if check_new_error_flag% an error flag needs to be created in the current cell


 93  ProjData.(FFName)= ~Mask;


 94  ProjData.VarDimName=[FieldData.VarDimName FieldData.VarDimName(CellInfo{icell}.VarIndex(1))];


 95  else


[672]  96  ProjData.(FFName)=FieldData.(FFName)  ~Mask;


 97  end


[639]  98  end


 99  end


 100  end


 101 


[672]  102 

