source: issm/oecreview/Archive/14064-14311/ISSM-14194-14195.diff

Last change on this file was 14312, checked in by Mathieu Morlighem, 12 years ago

Added 14064-14311

File size: 6.4 KB
  • ../trunk-jpl/src/m/plot/plot_gridded.m

     
    2323%Interpolating data on grid
    2424[x_m y_m data_grid]=InterpFromMeshToGrid(elements,x,y,data,xlim(1),ylim(2),post,post,round(diff(ylim)/post),round(diff(xlim)/post),NaN);
    2525if size(data_grid,1)<3 | size(data_grid,2)<3,
    26         error('data_grid size too small in plot_gridded, check posting and uni');
     26        error('data_grid size too small in plot_gridded, check posting and units');
    2727end
    2828
    2929%Get and change colormap
  • ../trunk-jpl/src/m/plot/kmlgridded.m

     
     1function kmlgridded(md,data,varargin)
     2
     3%process options
     4options=pairoptions(varargin{:});
     5
     6%process options
     7options=changefieldvalue(options,'coord','latlon');
     8
     9%process mesh and data
     10[x y z elements is2d isplanet]=processmesh(md,[],options);
     11[data datatype]=processdata(md,data,options);
     12
     13%check is2d
     14if ~is2d,
     15        error('buildgridded error message: gridded not supported for 3d meshes, project on a layer');
     16end
     17
     18%Get xlim and ylim (used to extract radar image)
     19xlim=[min(x) max(x)];
     20ylim=[min(y) max(y)];
     21post=getfieldvalue(options,'posting',diff(xlim)/1000);
     22if(diff(xlim)/post>10000),
     23        error(['posting too large']);
     24end
     25
     26%Interpolating data on grid
     27[x_m y_m data_grid]=InterpFromMeshToGrid(elements,x,y,data,xlim(1),ylim(2),post,post,round(diff(ylim)/post),round(diff(xlim)/post),NaN);
     28if size(data_grid,1)<3 | size(data_grid,2)<3,
     29        error('data_grid size too small, check posting and units');
     30end
     31pos=find(isinf(data_grid));
     32if ~isempty(pos),
     33        disp('Warning: removing Infs from vector (probably log(0)?)');
     34        data_grid(pos)=NaN;
     35end
     36
     37%Process data_grid: add white in NaN and correct caxis accordingly
     38data_nan=find(isnan(data_grid));
     39data_min=min(data_grid(:));
     40data_max=max(data_grid(:));
     41if exist(options,'caxis'),
     42        caxis_opt=getfieldvalue(options,'caxis');
     43        data_grid(find(data_grid<caxis_opt(1)))=caxis_opt(1);
     44        data_grid(find(data_grid>caxis_opt(2)))=caxis_opt(2);
     45        data_min=caxis_opt(1);
     46        data_max=caxis_opt(2);
     47end
     48
     49%Get colormap
     50colorm = getcolormap(options);
     51len    = size(colorm,1);
     52ind = ceil((len-1)*(data_grid-data_min)/(data_max - data_min + eps) +1);
     53ind(find(ind>len))=len;
     54ind(find(ind<1)  )=1;
     55ind(find(isnan(ind)))=1;
     56image_rgb=zeros(size(data_grid,1),size(data_grid,2),3);
     57r=colorm(:,1); image_rgb(:,:,1)=r(ind); clear r;
     58g=colorm(:,2); image_rgb(:,:,2)=g(ind); clear g;
     59b=colorm(:,3); image_rgb(:,:,3)=b(ind); clear b;
     60
     61%Deal with alpha
     62alpha=getfieldvalue(options,'alpha',.8);
     63alphaMatrix = alpha*ones(size(data_grid));
     64alphaMatrix(data_nan) = 0;
     65
     66%write kml
     67kmlfilename=getfieldvalue(options,'kmlfilename','tempfile.kml');
     68kmlroot=getfieldvalue(options,'kmlroot','./');
     69kmlimagename=getfieldvalue(options,'kmlimagename','tempimage');
     70kmlresolution=getfieldvalue(options,'kmlresolution',1);
     71kmlfolder=getfieldvalue(options,'kmlfolder','Ground Overlay');
     72kmlfolderdescription=getfieldvalue(options,'kmlfolderdescription','');
     73kmlgroundoverlayname=getfieldvalue(options,'kmlgroundoverlayname','ground overlay');
     74kmlgroundoverlaydescription=getfieldvalue(options,'kmlgroundoverlaydescription','description');
     75
     76%write png
     77imwrite(image_rgb,[kmlimagename '.png'],'png','alpha',alphaMatrix);
     78clear image_rgb alphaMatrix
     79
     80%prepare colorbar
     81iscolorbar=0;
     82if strcmpi(getfieldvalue(options,'colorbar','on'),'on'),
     83        X = linspace(0,1,len)';
     84        Xlab = round(linspace(data_min,data_max,len+1));
     85        html = ['<TABLE border=' num2str(1) ' bgcolor=#FFFFFF>',10];
     86
     87        for k=len:-1:1
     88                f = (Xlab(k)-data_min)/(data_max-data_min);
     89                if f<0, f=0; end
     90                if f>1, f=1; end
     91                polyColor(1,1) = interp1(X,colorm(:,1),f);
     92                polyColor(1,2) = interp1(X,colorm(:,2),f);
     93                polyColor(1,3) = interp1(X,colorm(:,3),f);
     94                polyColorStr(1:2) = dec2hex(round(polyColor(1)*255),2);
     95                polyColorStr(3:4) = dec2hex(round(polyColor(2)*255),2);
     96                polyColorStr(5:6) = dec2hex(round(polyColor(3)*255),2);
     97                html = [html,'<TR><TD width="15px" bgcolor=#',polyColorStr, '>&nbsp;</TD>','<TD bgcolor=#FFFFFF>'];
     98                if k==1
     99                        html=[html,'&lt;= ',num2str(Xlab(k),'%g')];
     100                elseif k==len
     101                        html=[html,'&gt;= ',num2str(Xlab(k),'%g')];
     102                else
     103                        html=[html,num2str(Xlab(k),'%g'),' to ',num2str(Xlab(k+1),'%g'),'</TD>'];
     104                end
     105                html = [html,'</TR>',10];
     106        end
     107        html = [html,'</TABLE>'];
     108        iscolorbar = 1;
     109end
     110
     111%now write kml file
     112fid=fopen([kmlroot '/' kmlfilename],'w');
     113fprintf(fid,'%s\n','<?xml version="1.0" encoding="UTF-8"?>');
     114fprintf(fid,'%s\n','<kml xmlns="http://earth.google.com/kml/2.1">');
     115fprintf(fid,'%s\n','<Document>');
     116fprintf(fid,'%s%s%s\n','<name>',kmlfilename,'</name>');
     117if iscolorbar,
     118        fprintf(fid,'<Placemark id="colorbar">\n');
     119        fprintf(fid,'%s%s%s\n','<name>','click the icon to see the colorbar','</name>');
     120        fprintf(fid,'%s%s%s\n','<description>','Ground overlay colorbar','</description>');
     121        fprintf(fid,'<visibility>1</visibility>\n');
     122        fprintf(fid,['<description>',10,'<![CDATA[' html ']]>',10,'</description>',10,'\n']);
     123        fprintf(fid,['<Style><IconStyle><scale>1</scale><Icon><href>http://maps.google.com/mapfiles/kml/shapes/donut.png</href></Icon></IconStyle><ListStyle></ListStyle></Style><Point id="poly_colorbar">\n']);
     124        fprintf(fid,'<altitudeMode>clampToGround</altitudeMode>\n');
     125        fprintf(fid,'<extrude>1</extrude>\n');
     126        fprintf(fid,'<tessellate>1</tessellate>\n');
     127        fprintf(fid,'%s%g,%g%s\n','<coordinates>',max(x),mean(y),'</coordinates>');
     128        fprintf(fid,'</Point>\n');
     129        fprintf(fid,'</Placemark>\n');
     130end
     131fprintf(fid,'%s\n','<GroundOverlay id="groundoverlay">');
     132fprintf(fid,'%s%s%s\n','<name>',kmlgroundoverlayname,'</name>');
     133fprintf(fid,'%s\n','<description>',kmlgroundoverlaydescription,'</description>');
     134fprintf(fid,'%s%s.%s%s\n','<Icon>',kmlimagename,'png','</Icon>');
     135fprintf(fid,'%s\n','<LatLonBox>');
     136fprintf(fid,'%s%f%s\n','<north>',max(y_m),'</north>');
     137fprintf(fid,'%s%f%s\n','<south>',min(y_m),'</south>');
     138fprintf(fid,'%s%f%s\n','<east>',max(x_m),'</east>');
     139fprintf(fid,'%s%f%s\n','<west>',min(x_m),'</west>');
     140fprintf(fid,'%s\n','<rotation>0</rotation>');
     141fprintf(fid,'%s\n','</LatLonBox>');
     142fprintf(fid,'%s\n','</GroundOverlay>');
     143fprintf(fid,'%s\n','</Document>');
     144fprintf(fid,'%s\n','</kml>');
     145fclose(fid);
Note: See TracBrowser for help on using the repository browser.