Матрица для скользящего окна

[Иван Стрельников]

Здравствуйте!
Возникла необходимость применить фильтр с радиальной маской (функция focal{raster}). Для построения матрицы использую следующее:

Код: Выделить всё

Circ_mask<-function(x){
	#Функция для построения матрицы радиального фильтра 
	#с заданным радиусом
	Matrix<-matrix(NA,1+(2*x),1+(2*x))
	Pos<-seq(-x,+x)
	for (i in c(1:nrow(Matrix))){
		for (j in c(1:ncol(Matrix))){
			#Расчет расстояния от центра к текущему 
			#положению
			VecL<-sqrt((Pos[i]**2)+(Pos[j]**2))
			if (x>=VecL){
				Matrix[i,j]<-1
			}
		}
	}
	return(Matrix)
}

Собственно вопрос -- есть ли в R готовая функция для подобных целей (сам пока не нашел)? Или, может, посоветуете более эффективный метод. Эта функция сильно грузная при больших радиусах. Плюс, если нужно будет добавить веса, то придется еще парится.
Ну и полностью тупой вопрос: если, например, я сделаю следующий вызов своей функции

Код: Выделить всё

> A<-Circ_mask(5)
> A
      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11]
 [1,]   NA   NA   NA   NA   NA    1   NA   NA   NA    NA    NA
 [2,]   NA   NA    1    1    1    1    1    1    1    NA    NA
 [3,]   NA    1    1    1    1    1    1    1    1     1    NA
 [4,]   NA    1    1    1    1    1    1    1    1     1    NA
 [5,]   NA    1    1    1    1    1    1    1    1     1    NA
 [6,]    1    1    1    1    1    1    1    1    1     1     1
 [7,]   NA    1    1    1    1    1    1    1    1     1    NA
 [8,]   NA    1    1    1    1    1    1    1    1     1    NA
 [9,]   NA    1    1    1    1    1    1    1    1     1    NA
[10,]   NA   NA    1    1    1    1    1    1    1    NA    NA
[11,]   NA   NA   NA   NA   NA    1   NA   NA   NA    NA    NA

, то какой у меня будет радиус скользящего окна 5 или, все же 5.5?

[gamm]

циклы в R писать вредно, там много чего есть вместо. Чтобы маску создать, первое, что на ум приходит - expand.grid()

Код: Выделить всё

p.R<-50
p.D<-2*p.R+1
p.grd<-expand.grid(-p.R:p.R,-p.R:p.R)
p.mask<-matrix(apply(p.grd,1,function(x) { if(sum(x^2) <= p.R^2) 1 else NA; }),ncol=p.D)

можно и быстрее, но и это сойдет. Веса очевидным образом вставляются в функцию.

Да, а радиус, естественно, 5, поскольку (5 - (-5))/2=5, а не 5.5

[KolesovDmitry]

Посмоторите пакет анализаизображений EBImage, там были фильтры, как предустановленные, так ипользовательские

[ Сообщение с мобильного устройства ]

[Иван Стрельников]

Большое спасибо, gamm!
Про for--циклы слышал несколько раз, но до сих пор, использовать Xapply в привычку не вошло -- буду исправляться. На счет радиуса, то в замешательство меня ввели размеры матрицы. Т.е. если рассматривать значения растра как значения в точках -- центрах пикселей, то все нормально. Но если воспринимать пиксели как четырехугольники, имеющие какую-то протяженность, то получается на 1 четырехугольник больше чем заданный диаметр.
По поводу EBImage, он вроде упразднен, но мануалы по нему в сети еще есть.

[gamm]

По поводу EBImage, он вроде упразднен, но мануалы по нему в сети еще есть.

он живет не в основном архиве CRAN, а в Bioconductor. И ориентирован на цитологию с соответствующими причудами.

Если уж чем-то и пользоваться для обработки изображений, то ITK, там вообще все есть. Для него есть пакет R - SimpleITK (живет на сайте ITK, нужно собирать самому), и есть другие проекты, в которые скорее R пристраивают к ITK, чем наоборот.

P.S. на растре все размеры считаем в ячейках.