R help - Oct 2010 - Contour Plot on a non Rectangular Grid

Dear All,
I would like to plot a scalar (e.g. a temperature) on a non-rectangular 
domain (or even better: I would simply like to be able to draw a contour 
plot on an arbitrary 2D domain). I wonder if there is any tool to 
achieve that with R. I did some online search in particular on the list 
archives, found several queries similar to this one but was not able to 
find any conclusive answer.
I am interested in the following 2 options

(1) just read a file of the form

x1  y1  z1
x2  y2  x2
... ... ...
xn  yn  zn

where the set of {xi} and {yi} are coordinates on an arbitrary domain 
and {zi} are the values of the scalar for the corresponding {x,y} 
coordinates.
(2) Sometimes the domain where I want to draw a contour plot is nothing 
too fancy and the scalar itself is given by an analytical function. 
Consider e.g. the case of a circle of radius R=pi/2 centered about the 
origin and a function like

z=f(x,y)=abs(cos(y))

NB: in this case a satisfactory solution could be to plot z on a 
rectangular grid and then clip a circular region
To fix the ideas, the final result in this case (with a colorjet map) 
should look like this

http://dl.dropbox.com/u/5685598/scalar_plot.pdf


Any suggestion is appreciated.
Many thanks

Lorenzo

On 24-Oct-10 11:30:57, Lorenzo Isella wrote:
> Dear All,
> I would like to plot a scalar (e.g. a temperature) on a non-rectangular
> domain (or even better: I would simply like to be able to draw a
> contour plot on an arbitrary 2D domain). I wonder if there is any
> tool to achieve that with R. I did some online search in particular
> on the list archives, found several queries similar to this one but
> was not able to find any conclusive answer.
> I am interested in the following 2 options
> 
> (1) just read a file of the form
> 
> x1  y1  z1
> x2  y2  x2
> ... ... ...
> xn  yn  zn
> 
> where the set of {xi} and {yi} are coordinates on an arbitrary domain 
> and {zi} are the values of the scalar for the corresponding {x,y} 
> coordinates.
> (2) Sometimes the domain where I want to draw a contour plot is nothing
> too fancy and the scalar itself is given by an analytical function. 
> Consider e.g. the case of a circle of radius R=pi/2 centered about the 
> origin and a function like
> 
> z=f(x,y)=abs(cos(y))
> 
> NB: in this case a satisfactory solution could be to plot z on a 
> rectangular grid and then clip a circular region
> To fix the ideas, the final result in this case (with a colorjet map) 
> should look like this
> 
> http://dl.dropbox.com/u/5685598/scalar_plot.pdf
> 
> Any suggestion is appreciated.
> Many thanks
> 
> Lorenzo
For your option (1), the fundamental issue is interpolation.
There are many methods for this, with different proprties!
An R Site Search on "interpolation" yields a lot of hits.
One (which is fairly basic, but may suit your purposes) is
the interpp() function in package akima:

http://finzi.psych.upenn.edu/R/library/akima/html/interpp.html

Hoping this helps,
Ted.

--------------------------------------------------------------------
E-Mail: (Ted Harding) <ted.harding at wlandres.net>
Fax-to-email: +44 (0)870 094 0861
Date: 24-Oct-10                                       Time: 12:51:03
------------------------------ XFMail ------------------------------

On 10/24/2010 01:51 PM, (Ted Harding) wrote:
> On 24-Oct-10 11:30:57, Lorenzo Isella wrote:
>> Dear All,
>> I would like to plot a scalar (e.g. a temperature) on a non-rectangular
>> domain (or even better: I would simply like to be able to draw a
>> contour plot on an arbitrary 2D domain). I wonder if there is any
>> tool to achieve that with R. I did some online search in particular
>> on the list archives, found several queries similar to this one but
>> was not able to find any conclusive answer.
>> I am interested in the following 2 options
>>
>> (1) just read a file of the form
>>
>> x1  y1  z1
>> x2  y2  x2
>> ... ... ...
>> xn  yn  zn
>>
>> where the set of {xi} and {yi} are coordinates on an arbitrary domain
>> and {zi} are the values of the scalar for the corresponding {x,y}
>> coordinates.
>> (2) Sometimes the domain where I want to draw a contour plot is nothing
>> too fancy and the scalar itself is given by an analytical function.
>> Consider e.g. the case of a circle of radius R=pi/2 centered about the
>> origin and a function like
>>
>> z=f(x,y)=abs(cos(y))
>>
>> NB: in this case a satisfactory solution could be to plot z on a
>> rectangular grid and then clip a circular region
>> To fix the ideas, the final result in this case (with a colorjet map)
>> should look like this
>>
>> http://dl.dropbox.com/u/5685598/scalar_plot.pdf
>>
>> Any suggestion is appreciated.
>> Many thanks
>>
>> Lorenzo
>
> For your option (1), the fundamental issue is interpolation.
> There are many methods for this, with different proprties!
> An R Site Search on "interpolation" yields a lot of hits.
> One (which is fairly basic, but may suit your purposes) is
> the interpp() function in package akima:
>
> http://finzi.psych.upenn.edu/R/library/akima/html/interpp.html
>
> Hoping this helps,
> Ted.
>
> --------------------------------------------------------------------
> E-Mail: (Ted Harding)<ted.harding at wlandres.net>
> Fax-to-email: +44 (0)870 094 0861
> Date: 24-Oct-10                                       Time: 12:51:03
> ------------------------------ XFMail ------------------------------
Hi,
And thanks for helping. I am anyway a bit puzzled, since case (1) is not 
only a matter of interpolation. Probably the point I did not make clear 
(my fault) is that case (1) in my original email does not refer to an 
irregular grid on a rectangular domain; the set of (x,y) coordinate 
could stand e.g. a flat metal slab along which I have temperature 
measurements. The slab could be e.g. elliptical or any other funny 
shape. What also matters is that the final outcome should not look 
rectangular, but by eye one should be able to tell the shape of the slab.
Case (1) is a generalization of case (2) where I do not have either an 
analytical expression for the surface not for the scalar.
Cheers

Lorenzo

On 24.10.2010 14:14, Lorenzo Isella wrote:
> On 10/24/2010 01:51 PM, (Ted Harding) wrote:
>> On 24-Oct-10 11:30:57, Lorenzo Isella wrote:
>>> Dear All,
>>> I would like to plot a scalar (e.g. a temperature) on a
non-rectangular
>>> domain (or even better: I would simply like to be able to draw a
>>> contour plot on an arbitrary 2D domain). I wonder if there is any
>>> tool to achieve that with R. I did some online search in particular
>>> on the list archives, found several queries similar to this one but
>>> was not able to find any conclusive answer.
>>> I am interested in the following 2 options
>>>
>>> (1) just read a file of the form
>>>
>>> x1 y1 z1
>>> x2 y2 x2
>>> ... ... ...
>>> xn yn zn
>>>
>>> where the set of {xi} and {yi} are coordinates on an arbitrary
domain
>>> and {zi} are the values of the scalar for the corresponding {x,y}
>>> coordinates.
>>> (2) Sometimes the domain where I want to draw a contour plot is
nothing
>>> too fancy and the scalar itself is given by an analytical function.
>>> Consider e.g. the case of a circle of radius R=pi/2 centered about
the
>>> origin and a function like
>>>
>>> z=f(x,y)=abs(cos(y))
>>>
>>> NB: in this case a satisfactory solution could be to plot z on a
>>> rectangular grid and then clip a circular region
>>> To fix the ideas, the final result in this case (with a colorjet
map)
>>> should look like this
>>>
>>> http://dl.dropbox.com/u/5685598/scalar_plot.pdf
>>>
>>> Any suggestion is appreciated.
>>> Many thanks
>>>
>>> Lorenzo
>>
>> For your option (1), the fundamental issue is interpolation.
>> There are many methods for this, with different proprties!
>> An R Site Search on "interpolation" yields a lot of hits.
>> One (which is fairly basic, but may suit your purposes) is
>> the interpp() function in package akima:
>>
>> http://finzi.psych.upenn.edu/R/library/akima/html/interpp.html
>>
>> Hoping this helps,
>> Ted.
>>
>> --------------------------------------------------------------------
>> E-Mail: (Ted Harding)<ted.harding at wlandres.net>
>> Fax-to-email: +44 (0)870 094 0861
>> Date: 24-Oct-10 Time: 12:51:03
>> ------------------------------ XFMail ------------------------------
>
> Hi,
> And thanks for helping. I am anyway a bit puzzled, since case (1) is not
> only a matter of interpolation. Probably the point I did not make clear
> (my fault) is that case (1) in my original email does not refer to an
> irregular grid on a rectangular domain; the set of (x,y) coordinate
> could stand e.g. a flat metal slab along which I have temperature
> measurements. The slab could be e.g. elliptical or any other funny
> shape. What also matters is that the final outcome should not look
> rectangular, but by eye one should be able to tell the shape of the slab.
> Case (1) is a generalization of case (2) where I do not have either an
> analytical expression for the surface not for the scalar.
> Cheers

What about the facilities in package rgl then?

Uwe Ligges

> Lorenzo
>
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.

On Oct 24, 2010, at 4:30 AM, Lorenzo Isella wrote:
> Dear All,
> I would like to plot a scalar (e.g. a temperature) on a non- 
> rectangular domain (or even better: I would simply like to be able  
> to draw a contour plot on an arbitrary 2D domain). I wonder if there  
> is any tool to achieve that with R. I did some online search in  
> particular on the list archives, found several queries similar to  
> this one but was not able to find any conclusive answer.
One implemented approach to this exists with the rms/Hmisc package  
combination. The perimeter function is used to define a region within  
which the are a sufficient number of cases and the perimeter object is  
passed to the bplot function, which is a wrapper for a lattice  
contourplot call. There is no reason you couldn't emulate

> I am interested in the following 2 options
>
> (1) just read a file of the form
>
> x1  y1  z1
> x2  y2  x2
> ... ... ...
> xn  yn  zn
>
> where the set of {xi} and {yi} are coordinates on an arbitrary  
> domain and {zi} are the values of the scalar for the corresponding  
> {x,y} coordinates.
> (2) Sometimes the domain where I want to draw a contour plot is  
> nothing too fancy and the scalar itself is given by an analytical  
> function. Consider e.g. the case of a circle of radius R=pi/2  
> centered about the origin and a function like
>
> z=f(x,y)=abs(cos(y))
>
That defines the contours but does not restrict the domain.
> NB: in this case a satisfactory solution could be to plot z on a  
> rectangular grid and then clip a circular region
> To fix the ideas, the final result in this case (with a colorjet  
> map) should look like this
>
> http://dl.dropbox.com/u/5685598/scalar_plot.pdf
And that color encoded output would not be the output of a contourplot  
but is more like a levelplot or an image plot. Nonetheless, the  
perimeter and bplot combination can deliver a similar result if you  
supply either code or data as a suitable test case for analysis and  
display.

>> Hi,
>> And thanks for helping. I am anyway a bit puzzled, since case (1) is
not
>> only a matter of interpolation. Probably the point I did not make clear
>> (my fault) is that case (1) in my original email does not refer to an
>> irregular grid on a rectangular domain; the set of (x,y) coordinate
>> could stand e.g. a flat metal slab along which I have temperature
>> measurements. The slab could be e.g. elliptical or any other funny
>> shape. What also matters is that the final outcome should not look
>> rectangular, but by eye one should be able to tell the shape of the
slab.
>> Case (1) is a generalization of case (2) where I do not have either an
>> analytical expression for the surface not for the scalar.
>> Cheers
>
>
> What about the facilities in package rgl then?
>
> Uwe Ligges
Hello,
I feel I am drowning in a glass of water.
Consider the following snippet at the end of the email, where I 
generated a set of {x,y,s=f(x,y)} values, i.e. a set of 2D coordinates + 
a scalar on a circle.
Now, I can get a scatterplot in 3D, but how to get a 2D surface 
plot/levelplot?
An idea could be to artificially set the z coordinate of the plot as a 
constant (instead of having it equal to s as in the scatterplot) and 
calculate the colormap with the values of s, along the lines of the 
volcano example + surface plot at

http://bit.ly/9MRncd

but I am experiencing problems. However, should I really go through all 
this? There is nothing truly 3D in the plot that I have in mind, you can 
think of it as e.g. some temperature measurement along a tube cross section.
Any help is appreciated.
Cheers

Lorenzo



########################################
library(scatterplot3d)
library(rgl)


R <- pi/2

n <- 100

x <- y <- seq(-R,R, length=n)

xys <- c()

temp <- seq(3)

for (i in seq(n)){

for (j in seq(n))

#check I am inside the circle
   if ((sqrt(x[i]^2+y[j]^2))<=R){

temp[1] <- x[i]
temp[2] <- y[j]
temp[3] <- abs(cos(y[j]))
xys <- rbind(xys,temp)

   }


}




scatterplot3d(xys[,1], xys[,2], xys[,3], highlight.3d=TRUE,
  col.axis="blue", col.grid="lightblue",
  main="scatterplot3d - 2", pch=20)



#########################################################

On Oct 24, 2010, at 9:30 AM, Lorenzo Isella wrote:
>
>>> Hi,
>>> And thanks for helping. I am anyway a bit puzzled, since case (1)  
>>> is not
>>> only a matter of interpolation. Probably the point I did not make  
>>> clear
>>> (my fault) is that case (1) in my original email does not refer to
>>> an
>>> irregular grid on a rectangular domain; the set of (x,y) coordinate
>>> could stand e.g. a flat metal slab along which I have temperature
>>> measurements. The slab could be e.g. elliptical or any other funny
>>> shape. What also matters is that the final outcome should not look
>>> rectangular, but by eye one should be able to tell the shape of  
>>> the slab.
>>> Case (1) is a generalization of case (2) where I do not have  
>>> either an
>>> analytical expression for the surface not for the scalar.
>>> Cheers
>>
>>
>> What about the facilities in package rgl then?
>>
>> Uwe Ligges
>
> Hello,
> I feel I am drowning in a glass of water.
Not sure what we are supposed to make of this.
> Consider the following snippet at the end of the email, where I  
> generated a set of {x,y,s=f(x,y)} values, i.e. a set of 2D  
> coordinates + a scalar on a circle.
> Now, I can get a scatterplot in 3D, but how to get a 2D surface plot/ 
> levelplot?
You were advised to look at rms. Why have you dismissed this  
suggestion? Using your data setup below and packaging into a dataframe.

require(rms)
  ddf <- datadist(xysf <- as.data.frame(xys))
  olsfit <- ols(V3~rcs(V1,3)+rcs(V2,3), data=xysf)

  bounds <- perimeter(xysf$V1, xysf$V2)
  plot(xysf$V1, xysf$V2) #demonstrates the extent of the data
  bplot(Predict(olsfit, V1,V2),  perim=bounds)  # a levelplot is the  
default

  bplot(Predict(olsfit, V1,V2),  perim=bounds, lfun=contourplot)
  bplot(Predict(olsfit, V1,V2),  perim=bounds, lfun=contourplot,  
xlim=c(-2.5,2.5))
# to demonstrate that perimeter works

# and as expected this shows very little variability d/t V1
  olsfit  # note that
 > anova(olsfit)
                 Analysis of Variance          Response: V3

  Factor          d.f. Partial SS   MS           F          P
  V1                 2   0.01618738 8.093691e-03      19.47 <.0001
   Nonlinear         1   0.01618738 1.618738e-02      38.93 <.0001
  V2                 2 470.67057254 2.353353e+02  566040.95 <.0001
   Nonlinear         1 470.67057254 4.706706e+02 1132081.91 <.0001
  TOTAL NONLINEAR    2 527.78127558 2.638906e+02  634723.80 <.0001
  REGRESSION         4 527.78127558 1.319453e+02  317361.90 <.0001
  ERROR           7663   3.18594315 4.157566e-04
# most the the regression SS is in the V2 variable
# Q.E.D.
-- 
David,
> An idea could be to artificially set the z coordinate of the plot as  
> a constant (instead of having it equal to s as in the scatterplot)  
> and calculate the colormap with the values of s, along the lines of  
> the volcano example + surface plot at
>
> http://bit.ly/9MRncd
>
> but I am experiencing problems. However, should I really go through  
> all this? There is nothing truly 3D in the plot that I have in mind,  
> you can think of it as e.g. some temperature measurement along a  
> tube cross section.
> Any help is appreciated.
> Cheers
>
> Lorenzo
>
> ########################################
> library(scatterplot3d)
> library(rgl)
> R <- pi/2
>
> n <- 100
>
> x <- y <- seq(-R,R, length=n)
>
> xys <- c()
>
> temp <- seq(3)
>
> for (i in seq(n)){
>
>     for (j in seq(n))
>
> #check I am inside the circle
>        if ((sqrt(x[i]^2+y[j]^2))<=R){
>
>       temp[1] <- x[i]
>       temp[2] <- y[j]
>       temp[3] <- abs(cos(y[j]))
>      xys <- rbind(xys,temp)
>  }
> }
> scatterplot3d(xys[,1], xys[,2], xys[,3], highlight.3d=TRUE,
> col.axis="blue", col.grid="lightblue",
> main="scatterplot3d - 2", pch=20)
>
>
>
> #########################################################