# Reading JG Kunkel Eclogite Data Set 4.  JG Kunkel
smoo9<- function(M) {
	r<- length(M[,1]); c<-length(M[1,])
	rm1<- r-1; cm1<- c-1; rm2<- r-2; cm2<- c-2
	Mavg<- (8*M[2:rm1, 2:cm1] + 
	       M[1:rm2, 1:cm2] + 
	       M[1:rm2, 2:cm1] + 
	       M[1:rm2, 3:c] + 
	       M[2:rm1, 1:cm2] + 
	       M[2:rm1, 3:c] + 
	       M[3:r, 1:cm2] + 
	       M[3:r, 2:cm1] + 
	       M[3:r, 3:c])/16
	       M[2:rm1, 2:cm1]<-Mavg
	       Mavg<-M
	       Mavg
	                 }

filnam= "JK depression 1"
filnamCa<- paste(filnam,"_2.txt", sep='')
filnamP<- paste(filnam,"_3.txt", sep='')
filnamP2<- paste(filnam,"_1.txt", sep='')  # F-
filnamCl<- paste(filnam,"_4.txt", sep='')
XX<-read.table(filnamP, skip=25, nrows=2)
Ncols<- XX[1,4]
Nrows<- XX[2,4]
iNcols<- Ncols:1
iNrows<- Nrows:1
rm(XX)
quartz()
XXp<- as.matrix(read.table(filnamP, skip=33, nrows=Nrows, sep='\t'))
image(XXp, col=gray(0:255/255), main="P")
quartz()
XXp2<- as.matrix(read.table(filnamP2, skip=33, nrows=Nrows, sep='\t'))
image(XXp2, col=gray(0:255/255), main="F")
quartz()
# XXp<- (XXp + XXp2)/2
image(XXp2, col=gray(0:255/255), main="Pavg")
quartz()
XXCa<- as.matrix(read.table(filnamCa, skip=33, nrows=Nrows, sep='\t'))
image(XXCa, col=gray(0:255/255), main="Ca")
quartz()
XXCl<- as.matrix(read.table(filnamCl, skip=33, nrows=Nrows, sep='\t'))
image(XXCl, col=gray(0:255/255), main="Cl")
quartz()
XXCl<- smoo9(XXCl)
image(XXCl, col=gray(0:255/255), main="Cl smoothed")

require(pixmap)
XX3<- 0*XXp
quartz()
rat<- 40      
rat2<- 40
z <- pixmapRGB(c(rat*XXCl,XXCa,rat2*XXp), nrow=Nrows, ncol=Ncols, bbox=c(-1,-1,1,1))
plot(z, axes=TRUE, main="Red = Chloride; Green = calcite; Violet = tetraCaP")
write.pnm(z, file=paste(filnam,"_RGB.pnm", sep=''))