Best way to plot a surface whose matrix has values of +/- infinity

Hi there! What is the best way to plot a 3D surface whose matrix C = A ./ B contains values of +/- infinity (whenever entries in B are zero)? Currently I am using surf( ), which results in weird vertical cross sections. I looked into the fillMissing( ) and scatteredInterpolant( ) functions. Is there a best / recommended way to proceed, to plot a smooth surface? Thanks in advance,

[Torsten]

Do you mean C = A/B or C = A./B ?

[Noob]

Hi Torsten! ./

Answer 1

[Star Strider]

Without the actual matrices (or representative matrices), it is diifficult to determine what is best.

Note that fillmissing will only interpolate NaN values, not ±Inf since those aren’t actually ‘missing’. If you want to fill the ±Inf values to a specific value or simply interpolate them, first change the ±Inf values to NaN, and then use fillmissing with the appropriatee options, depending on what you want to fill the values with.

C = rand(10);
C(randi(numel(C),5,1)) = Inf;
disp(C)

    0.9588    0.5339    0.0804    0.2514    0.9898    0.6714    0.1631    0.9184    0.4323    0.1071
    0.1392       Inf    0.7912    0.7563    0.1575    0.9422    0.5456    0.6430       Inf    0.7044
    0.3319    0.2681    0.3178    0.2085    0.2741    0.3512    0.3111    0.8425    0.3927    0.8152
       Inf    0.6310    0.1859       Inf    0.8198    0.6044    0.2868    0.6497    0.2099    0.6087
    0.4576    0.6661    0.5924    0.4996    0.2698    0.0185    0.2673    0.7120    0.4400    0.4482
    0.9765    0.9727    0.3978    0.5604    0.9376    0.0350    0.3679    0.4305    0.6991    0.0103
    0.9861    0.9224    0.5639    0.3243    0.6420    0.0775    0.5452    0.5099    0.2087    0.8130
    0.2795    0.8829    0.9604    0.6973    0.7640    0.9333    0.2153    0.5529    0.6946    0.6082
    0.9332       Inf    0.2878    0.9074    0.2452    0.4789    0.6202    0.8002    0.1871    0.7009
    0.1150    0.5625    0.2731    0.5491    0.8665    0.0998    0.9425    0.4497    0.5860    0.4568

C(~isfinite(C)) = NaN

C = 10×10

0.9588 0.5339 0.0804 0.2514 0.9898 0.6714 0.1631 0.9184 0.4323 0.1071 0.1392 NaN 0.7912 0.7563 0.1575 0.9422 0.5456 0.6430 NaN 0.7044 0.3319 0.2681 0.3178 0.2085 0.2741 0.3512 0.3111 0.8425 0.3927 0.8152 NaN 0.6310 0.1859 NaN 0.8198 0.6044 0.2868 0.6497 0.2099 0.6087 0.4576 0.6661 0.5924 0.4996 0.2698 0.0185 0.2673 0.7120 0.4400 0.4482 0.9765 0.9727 0.3978 0.5604 0.9376 0.0350 0.3679 0.4305 0.6991 0.0103 0.9861 0.9224 0.5639 0.3243 0.6420 0.0775 0.5452 0.5099 0.2087 0.8130 0.2795 0.8829 0.9604 0.6973 0.7640 0.9333 0.2153 0.5529 0.6946 0.6082 0.9332 NaN 0.2878 0.9074 0.2452 0.4789 0.6202 0.8002 0.1871 0.7009 0.1150 0.5625 0.2731 0.5491 0.8665 0.0998 0.9425 0.4497 0.5860 0.4568

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

●

C = fillmissing(C, 'linear')

C = 10×10

0.9588 0.5339 0.0804 0.2514 0.9898 0.6714 0.1631 0.9184 0.4323 0.1071 0.1392 0.4010 0.7912 0.7563 0.1575 0.9422 0.5456 0.6430 0.4125 0.7044 0.3319 0.2681 0.3178 0.2085 0.2741 0.3512 0.3111 0.8425 0.3927 0.8152 0.3948 0.6310 0.1859 0.3540 0.8198 0.6044 0.2868 0.6497 0.2099 0.6087 0.4576 0.6661 0.5924 0.4996 0.2698 0.0185 0.2673 0.7120 0.4400 0.4482 0.9765 0.9727 0.3978 0.5604 0.9376 0.0350 0.3679 0.4305 0.6991 0.0103 0.9861 0.9224 0.5639 0.3243 0.6420 0.0775 0.5452 0.5099 0.2087 0.8130 0.2795 0.8829 0.9604 0.6973 0.7640 0.9333 0.2153 0.5529 0.6946 0.6082 0.9332 0.7227 0.2878 0.9074 0.2452 0.4789 0.6202 0.8002 0.1871 0.7009 0.1150 0.5625 0.2731 0.5491 0.8665 0.0998 0.9425 0.4497 0.5860 0.4568

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

Another option is to do someething similar with ‘B’ except to change the zero values to NaN and then use fillmissinng to interpolate them. Then do the division.

.

[Noob]

Hi Star Strider!

Does fillmissing( ) work better at interpolation if I increase the size of the matrices, i.e. refine the meshgrid? I have little-to-no numerical analysis background but I would guess that it helps. But for now, let me proceed with my current workflow to see how things look initially. Thanks!

[Star Strider]

My pleasure!

The fillmissing function works the same way regardless of the matriix size.

Another option is fillmissing2 (introduced in R2023a) that might provide better results.

C = randn(10);
C(randi(numel(C),25,1)) = Inf;
disp(C)

       Inf       Inf    0.3098   -0.0164       Inf   -1.1556    1.6242   -0.1819   -1.0374   -0.5292
    0.2479   -1.1064    1.2551       Inf   -0.2629    0.8479    0.2331       Inf   -0.8725    0.5190
       Inf    1.6133   -1.2166   -0.1762   -0.6104   -0.6279    0.0181   -0.1940       Inf    0.0363
       Inf    0.6911       Inf    0.4206    0.0051    0.3790    0.9779    1.0189   -0.3138   -0.7988
    0.2660   -3.2065   -0.6493   -1.4479   -1.0254       Inf    0.6221   -0.5161    0.7479       Inf
   -0.1340    0.5418    1.0401       Inf       Inf   -0.0150   -1.2801       Inf    1.5681       Inf
   -1.5790   -0.0629    0.2349    0.1572   -0.6514    2.0185    1.3684    1.0410   -0.7846       Inf
    1.6414   -0.6006       Inf    0.1635   -1.0758   -1.1217   -0.4134    0.2934    0.3222    0.6822
       Inf       Inf    0.0813    0.2288   -0.6510    0.4778    0.2747       Inf       Inf   -0.6654
    0.5429    0.0018   -2.8620   -0.3351   -0.5000       Inf    1.2599       Inf   -1.2694   -0.5554

missingLocs = ~isfinite(C)

missingLocs = 10x10 logical array
   1   1   0   0   1   0   0   0   0   0
   0   0   0   1   0   0   0   1   0   0
   1   0   0   0   0   0   0   0   1   0
   1   0   1   0   0   0   0   0   0   0
   0   0   0   0   0   1   0   0   0   1
   0   0   0   1   1   0   0   1   0   1
   0   0   0   0   0   0   0   0   0   1
   0   0   1   0   0   0   0   0   0   0
   1   1   0   0   0   0   0   1   1   0
   0   0   0   0   0   1   0   1   0   0

figure
surfc(C)
colormap(turbo)
title('Original')

●

C1 = C;
C1(~isfinite(C1)) = NaN;
C1 = fillmissing(C1, 'nearest', MissingLocations=missingLocs)

C1 = 10×10

0.2479 -1.1064 0.3098 -0.0164 -0.2629 -1.1556 1.6242 -0.1819 -1.0374 -0.5292 0.2479 -1.1064 1.2551 -0.1762 -0.2629 0.8479 0.2331 -0.1940 -0.8725 0.5190 0.2479 1.6133 -1.2166 -0.1762 -0.6104 -0.6279 0.0181 -0.1940 -0.3138 0.0363 0.2660 0.6911 -0.6493 0.4206 0.0051 0.3790 0.9779 1.0189 -0.3138 -0.7988 0.2660 -3.2065 -0.6493 -1.4479 -1.0254 -0.0150 0.6221 -0.5161 0.7479 -0.7988 -0.1340 0.5418 1.0401 0.1572 -0.6514 -0.0150 -1.2801 1.0410 1.5681 0.6822 -1.5790 -0.0629 0.2349 0.1572 -0.6514 2.0185 1.3684 1.0410 -0.7846 0.6822 1.6414 -0.6006 0.0813 0.1635 -1.0758 -1.1217 -0.4134 0.2934 0.3222 0.6822 0.5429 0.0018 0.0813 0.2288 -0.6510 0.4778 0.2747 0.2934 -1.2694 -0.6654 0.5429 0.0018 -2.8620 -0.3351 -0.5000 0.4778 1.2599 0.2934 -1.2694 -0.5554

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

figure
surfc(C1)
colormap(turbo)
title('fillmissing')

C2 = C;
C2(~isfinite(C2)) = NaN;
C2 = fillmissing2(C2, 'nearest', MissingLocations=missingLocs)

C2 = 10×10

0.2479 0.3098 0.3098 -0.0164 -1.1556 -1.1556 1.6242 -0.1819 -1.0374 -0.5292 0.2479 -1.1064 1.2551 -0.2629 -0.2629 0.8479 0.2331 -0.8725 -0.8725 0.5190 1.6133 1.6133 -1.2166 -0.1762 -0.6104 -0.6279 0.0181 -0.1940 0.0363 0.0363 0.6911 0.6911 0.4206 0.4206 0.0051 0.3790 0.9779 1.0189 -0.3138 -0.7988 0.2660 -3.2065 -0.6493 -1.4479 -1.0254 0.6221 0.6221 -0.5161 0.7479 -0.7988 -0.1340 0.5418 1.0401 0.1572 -0.0150 -0.0150 -1.2801 1.5681 1.5681 1.5681 -1.5790 -0.0629 0.2349 0.1572 -0.6514 2.0185 1.3684 1.0410 -0.7846 0.6822 1.6414 -0.6006 0.1635 0.1635 -1.0758 -1.1217 -0.4134 0.2934 0.3222 0.6822 0.5429 0.0813 0.0813 0.2288 -0.6510 0.4778 0.2747 0.2934 -0.6654 -0.6654 0.5429 0.0018 -2.8620 -0.3351 -0.5000 1.2599 1.2599 -1.2694 -1.2694 -0.5554

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

figure
surfc(C2)
colormap(turbo)
title('fillmissing2')

After some experimentatiion, the 'nearest' method is probably the most reliable.

.

[Star Strider]

I don’t have your data, so I can’t determine what the problem is with it.

Perhaps this would work better:

B(abs(B) <= 1E-15) = NaN

It first creates a logical matrix with ‘(abs(B) <= 1E-15)’ and sets the corresponding true values of ‘B’ to NaN.

Try that.

.

[Noob]

Ok, I'll try that. I've only recently started using logical matrices. Let me proceed with the fillMissing2 function after. Thanks!

[Noob]

Hi Star Strider,

I get this warning, when using fillMissing2( ):

Warning: The underlying triangulation is empty - the points may be collinear.

[Star Strider]

Again, since I don¹t have your data, I can’t offer any specific suggestions on that, other than to use fillmissing instead. It should give decent results, since it did in my example matrices in my previous Comment. The fillmissing function won’t care about the spatial characteristics, since it operates on the linear indices, vector-wise, not spatially.

That may be the best that you can do.

[Noob]

Hi Star Strider,

The method "movmedian",3 seems to work well. No more errors about triangulation. Let me check things out. Thanks again!

[Star Strider]

My pleasure!

[Noob]

Hi Star Strider,

Many of the NaN values are filled in, but some NaN values are not filled in. I'll try to see what's going on. Matlab's not giving me an error message though. I also found that when using a bigger number, like "movmedian",10, more of the NaN entries get filled in.

[Noob]

fillmissing vs fillmssing2 is slightly interesting to see. I will post a new question: about a lingering vertical cross section in this surface that seems plainly wrong to me, as it would seem to indicate that the function is multi-valued, when it is not. Thanks again! Let me accept your answer now.

[Star Strider]

I’m not exactly certain what you’re doing at this point.

I discovered that as well, and solved it by telling fillmissing or fillmissing2 where the missing values are.

In my code, I did that with:

missingLocs = ~isfinite(C)

and then used that value (logical matrix) with MissingLocations to both function calls:

C1 = fillmissing(C1, 'nearest', MissingLocations=missingLocs)

and:

C2 = fillmissing2(C2, 'nearest', MissingLocations=missingLocs)

That worked.

.

[Star Strider]

As always, my pleasure!

[Noob]

Interesting! Let me try specifying MissingLocations now. Thanks!

[Star Strider]

As always, my pleasure!

[Noob]

Hi Star Strider!

I'm happy to report that I can now remove those vertical cross sections that were based on the function oscillating from -infinity to +infinity, due to function discontinuity. The problem now, having a detached surface, is to interpolate function values, to get a whole, connected, smooth surface, I think.

fillMissing( ) doesn't interpolate over all the NaN values in the Z matrix. So, as you said, I have to tell fillMissing( ) where it missed some NaNs, but I don't quite understand your code.

I will look at it again soon.

Thanks again for your help!

[Noob]

Hi Star Strider,

I just figured it out.

Thanks again!

[Star Strider]

As always, my pleasure!

[Noob]

Hi Star Strider,

Using a very crude bound gets me many more initial NaNs, and then the fully-interpolated surface looks elegant and simple, essentially matching a simple version of the math model I had in mind. Using a finer and finer bound gets me less initial NaNs, and a more complicated surface that, once fully interpolated by fillMissing( ), somewhat resembles the truer complicated surface associated with the more complicated math model.

Where should we draw the line when setting such a bound?

Should we aim to get the simplest figure, or to get the truest figure?

I'm not sure, but I've started to think in this direction.

Thanks!

[Walter Roberson]

In my opinion, if the function is oscillating to +/- inf, then it should be left that way, with no interpolation.

My opion might be a little different if the +/- inf values are caused by overflow of numeric representation, with the actual values being "reasonably" finite (if, for example, symbolic mathematics were used.)

Any theoretical +/- inf should be left as +/- inf, in my opinion. If that causes problems with graphing, then either "Tough luck" or else use ylim() / zlim() to cut off the infinite graph.

[Noob]

Hi @Walter Roberson! Yes, at the moment, I think I would agree with you. The spikes in function behavior should be left in as they are -- and explained in the main text / captions. I think for now, I will stop work on interpolation, etc. and get back to work on other things. Thanks for sharing!

[Star Strider]

@Noob — I’m not certain what you mean by ‘bound’. It seems to me that a more appropriate term might be ‘resolution’, in that a finer resolution (more data poiints within specific liimits) could be what you intend. In that regard, the best tradeoff between computation time and available memory might be best, and I cannot determine what that is because I have no idea what you are actually doing.

I would go for the most accurate depiction of the model, bearing in mind that the ±Inf values should probably still look like peaks or valleys, although with finite (flat) amplitudes with an appropriate value.

[Noob]

Hi Star Strider! Ok, got it; ultimately, I don't have the final say in such things. I can only propose ideas, and keep working, and moving forward. This particular figure and its associated model is undoubtedly difficult to understand. I want to comment on your answer about color shading, so I'll go there now -- thanks again!

[Star Strider]

As always, my pleasure!