Re: Problem with griddata and contourplot

[Ben Abbott]

On Dec 11, 2011, at 1:40 PM, Ben Abbott wrote:

> On Dec 9, 2011, at 11:53 AM, Walter White wrote:
> 
>> 2011/12/9 Ben Abbott <address@hidden>
>> 
>>> On Dec 9, 2011, at 11:04 AM, Walter White wrote:
>>> 
>>>> Hello,
>>>> 
>>>> I came across a problem and hope that you can help me.
>>>> 
>>>> I was trying to produce a contour plot of scattered x,y,z data,
>>>> therefore I was using griddata. So far everything works fine,
>>>> but when my mesh is too fine the following error message
>>>> appears and the script quits:
>>>> 
>>>> gnuplot (as of v4.2) only supports 3D filled triangular patches
>>>> 
>>>> Since the data is from a FEM software I have to plot it
>>>> using a fine mesh to see the details.
>>>> Does someone know a solution or a workaround to that?
>>>> 
>>>> Kind regards,
>>>> Walter
>>> 
>>> Can you tell us what commands you used?
>>> 
>>> And, if possible, provide a simple example for us to try ourselves?
>>> 
>>> Ben
>>> 
>>> <my_FEM_matrix.dat><FEM_matrix_stepwidth_1.00.png><FEM_matrix_stepwidth_0.50.png><FEM_matrix_stepwidth_0.30.png><FEM_matrix_stepwidth_0.10.png><FEM_matrix_stepwidth_0.05.png>
>> 
>> Thanks for the quick reply!
>> 
>> I used
>> 
>> ----------
>> FEM_matrix = load('./my_FEM_matrix.dat');
>> stepwidth = 0.03 %in m
>> x_density = min(FEM_matrix(:,1)) : stepwidth : max(FEM_matrix(:,1));
>> y_density = minFEM_matrix(:,2)) : stepwidth : max(FEM_matrix(:,2));
>> 
>> [XI,YI] = meshgrid(x_density,y_density);
>> 
>> x = FEM_matrix(:,1);
>> y = FEM_matrix(:,2);
>> z = FEM_matrix(:,3);
>> 
>> ZI = griddata(x,y,z,XI,YI);
>> 
>> contourf(XI,YI,ZI,50)
>> 
>> ylabel('y-axis [m]')
>> xlabel('x-axis [m])')
>> 
>> title('FEM matrix')
>> 
>> axis([min(FEM_matrix(:,1)) , max(FEM_matrix(:,1)) , min(FEM_matrix(:,2)) , 
>> max(FEM_matrix(:,2)) ]);
>> 
>> caxis([0 5])
>> hd = colorbar;
>> set(get(hd,'title'),'String','[ radian ]');
>> 
>> print('some_file_name.png', '-dpng','-S1024,768')
>> --------------------
>> 
>> I also tried to disable the call to axis(...), but this does not change a 
>> thing.
>> The function crashes with the error message at about stepwidth = 0.02.
>> 
>> ---------
>> error: gnuplot (as of v4.2) only supports 3D filled triangular patches
>> error: called from:
>> error:   /usr/share/octave/3.2.4/m/plot/__go_draw_axes__.m at line 567, 
>> column 6
>> error:   /usr/share/octave/3.2.4/m/plot/__go_draw_figure__.m at line 92, 
>> column 3
>> error:   /usr/share/octave/3.2.4/m/plot/gnuplot_drawnow.m at line 99, column 
>> 5
>> error:   /usr/share/octave/3.2.4/m/plot/print.m at line 707, column 7
>> -------
>> 
>> Another strange thing are the white areas on the plots that appear
>> is a finer mesh is used. I attached some pictures.
>> 
>> Kind regards,
>> Walter
> 
> There were some syntax problems with the script. In case anyone else would 
> like to take a look, I've corrected them (see below).
> 
> I've verified that Octave and Matlab produce the same ZI and that both 
> gnuplot and fltk produce similar results (different than matlab). The 
> comparative plots and the corrected script can be found at the link below.
> 
>       http://files.me.com/bpabbott/rd3bjp
> 
> The link will expire in 30 days.
> 
> I also compared the plots when contourf() is replaced by contour(). The 
> non-filled contour plots look essentially the same for gnuplot, fltk and 
> Matlab.
> 
> It appears that this is a bug in Octave's construction of the filled contour.
> 
> I've filed a bug report.
> 
>       https://savannah.gnu.org/bugs/index.php?35033
> 
> Ben

Walter,

As a work around, you can try imagesc(). Replace the call to contourf() with ...

        imagesc (XI, YI, ZI)

If you only want 50 levels, you can do something like ...

        Zmax = max (ZI(:));
        imagesc (XI, YI, Zmax*0.02*round(50*ZI/Zmax))

Ben