Pages

Showing posts with label Remote sensing. Show all posts
Showing posts with label Remote sensing. Show all posts

Saturday, 14 February 2015

MATLAB tutorial: Dividing image into blocks and applying a function

Often due to the limitation in memory, we want to divide an image into mxn blocks and process those blocks. If your block processing outputs an image then, you can use blockproc function in MATLAB. But, if processing function outputs points then one can’t use block processing. For that, one has to rely on image indexing to divide image into blocks. So here I will show you how you can divide image into blocks and process those individual blocks for any particular function using indexing. I am going to use block processing for finding pivot irrigation fields detection that was featured in this post. The trick is to apply matrix indexing to get chunk of block data and process them in a sequential manner. 

Here is the original image. We are going to divide it into mxn blocks. With the code you can specify any m or n value.

Original image

Block processing sequential numbering
With the following code, the block processing on individual blocks are performed. Feel free to copy the code and adapt it to your liking. Everything in the code is self-explanatory, so just go through the code line by line and you will understand whats going on.

The output of block processing for original image is as follows:

Block processing of detecting circles 

.

Wednesday, 14 January 2015

MATLAB Tutorial: Finding center pivot irrigation fields in a high resolution image

In this post, we will be talking about finding circular irrigation fields. The image was download from here. If you want to try, download it and play around. The approach I am using is detection of edges and finding edges that form circular edges. If you have heard of Hough’s transformation for detection of straight lines, the method is just an extension of it to detect circles. For Hough’s transformation, please go through here. If you have this excellent book, it has an elaborate explanation of Hough’s transformation to detect straight lines. Peoples have used Hough’s transformation for many different purpose. Many people use it for detecting straight building edges that can used to reconstruct buildings for 3-D buildings, generating city GML models etc.

For Hough’s transformation for straight line detection, you can use either PYTHON based Scikit-image or MATLAB. Circular Hough's transformation is also available in both scikit-image and openCV. The algorithm (imfindcircles) is available in MATLAB since 2013b version with image processing toolbox. I could not find the algorithm in any remote sensing software so far. My personnel view is that remote sensing software are very behind on incorporating state-of-the-art algorithms. So, knowing some coding either PYTHON, MATLAB or R will take you to greater heights in your professional path.

A circle is represented mathematically as:

(x-x_{center})^2 + (y - y_{center})^2 = r^2

where xcenter and ycenter are center of the circle and r is the radius of the circle. As you can see, there are 3 parameters to be fitted for cicular hough transformation.

Read the documentation of MATLAB, to get an idea about imfindcircles function parameters. So here basically, we are going to use imfindcircles function to detect center pivot irrigation fields in the image. In this image, there are only dark pivot circular fields surround by bright objects. So, we will be using only ‘dark’ mode of imfindcircles. The minimum circle and the maximum circle radii are image dependent so you need to provide those information with a little bit of data exploration.

Here is a code in MATLAB.

Original RGB image

RGB image with detected pivot irrigation fields
As you can see, 7 fields out of 9 fields were correctly detected. Two undetected fields( middle -top of the image) are also darker circles but due to its low contrast with surrounding, they were not detected. Even with the higher value for parameter 'sensitivity' and the low value of 'edgeThreshold' paramter, those two fields were undetected. You can further increse sensitivity and lower edgeThreshold parameters, to find those undetected circles but then you risk of finding many false alarms as well.

I have a gut feeling that with OBIA with ecognition, the process of finding center pivot irrigation fields would be not a straight forward procedure as with using the function  imfindcircles and the process would be much complex. Nevertheless, i will try it with eCognition in near future and report back.

Pixel-based based classification using any machine learning classifies will fail miserably for this case as center pivot irrigation fields ares spectrally similar to vegeation in other rectangular plots.

UPDATE


I spent some time in eCognition exploring a newer algorithm " template matching". The technique is not new but it has been incorporated with the last release of eCognition. The concept is given a template, the template moves over the image (a single layer) in a sliding window and calculates normalized cross-correlation simalarity between the template and the pixels within the sliding window. The result is an cross-correlation image. Subsequently a threshold value is used to find position of pixels with higher cross-correlation value.

The concept of normalized cross-correlation is shown below taken from a good presentation. Study the presentation in detail if you want.Notice border effect of the cross-correlation image below. This can be avoided if padding with replicated pixels are added in the image (commonly done in MATLAB) for any convolution procedure.

The detail explanation for performing the template matching in eCognition will follow sometime in future.

Cross-correlation concept

A genarated template with many samples.

Cross-Correlation image

Result of template matching in eCognition



Tuesday, 16 December 2014

eCognition Tutorial: Finding trees and buildings from LiDAR with limited information

I have worked a lot with LiDAR during my MSc time, but for my current work I am working more with optical images. I kind of miss LiDAR. When one of my friends came  to me with a problem related to LiDAR, I was very happy and decided to flaunt her my eCognition skills :-).

LiDAR data typically comes with many attributes like no. of returns, intensity, First Pulse Elevation, Last Pulse Elevation etc. Some data collectors even provides preliminary discrimination into ground and non-ground. But in this case, all we got is First Pulse data and Last Pulse data. And the desired output is discrimination between trees and buildings.

The data was borrowed from eCognition community. If someone has time to kill, head over there, get the data, roll up your sleeves and let’s do information extraction from LiDAR with eCognition, shall we?

My workflow:
  1. Create a diff image (FP-LP) and classify tree
    1. Assign pixels > 2 as high with multi-threshold segmentation
    2. Assign small objects < 10 pixels surrounded by high as high
    3. Perform opening and closing with ball shaped Structuring elements (SE). Opening step is required to remove footprint effect of LiDAR along the buildings edges.
    4. Assign objects with area > 20 pixels as tree
  2. Find buildings in Last Pulse images
    1. Perform chess-board segmentation with size 1 on unclassified objects
    2. Perform a Multi-Resolution Segmentation (MRS)
    3. Create a feature “Mean Difference to unclassified” feature within neighborhood of 20 pixels. For that a customized feature was created.
    4. Assign unclassified objects with Mean Difference to unclassified (20) > 4 m
    5. Merge buildings objects and assign small objects ( area < 100 pixels) classified as buildings as trees.
    6. Assign unclassified objects surrounded by buildings as buildings
The rule-set development took 20 minutes of my lunch time. The process takes 11 seconds for an area 360 m x 360 m and the result obtained is reasonably good. With a little more effort, the result can obviously be enhanced. Nevertheless, i gave myself a pat on the back.

When i will have more time in near future, i will compare the result with adopting a different methodology using lastools.

Rule-set in action

First Pulse
Diffrence image
Classification ( Yellow: Buildings, Green: Trees)

Last Pulse
UPDATE

Well, this morning my friend called me and told me that the classification of buildings is great. Can we get straight lines for buildings edges rather than zig-zag lines? Lets see what she is taking about. Yes, buildings are most of the times straight. But due to the effect of segmentation, our classified buildings edges are zig-zag.The problem can be tackled with native vector handling capability of eCognition. The algorithms were introduced in eCognition 9.0 that was released couple of months ago.

Zig-Zag edges problem of buildings
Approach 1

  • Convert building objects into a shp file
  • Use buiding orthogonalization algorithm (Chessboard: 7 pixels and Merge Threshold: 0.5)

As you can, the result is far from perfect.

Approach 2
  • Use mathematical morphology closing (SE: box 7x7 pixels) on building objects
  • Use mathematical morphology opening (SE: box 7x7 pixels) on building objects
  • Convert buildings objects into a shp file
  • Use buiding orthogonalization algorithm (Chessboard: 7 pixels and Merge Threshold: 0.5)

The result is much better than first approach. Here, the sequence of close-open must be followed since buildings objects that are loosely connected may break into separate objects if a sequence of open-close is adopted.

Boundary othrogonalization without Mathematical Morphology ( Yellow: building objects, Red: New boundary)
Boundary othrogonalization with Mathematical Morphology ( Yellow: building objects, Red: New boundary)




Tuesday, 18 November 2014

eCognition Tutorial: Find the closest classified object


Solution:

There are different ways to tackle this problem. First one should be aware that there are two different methods of computing image objects distance. Center of Gravity and Smallest Enclosing Rectangle. Scour the help doc to find what those two methods are and choose the method you want. You can set up this in the beginning of your rule set using rule set options algorithm. Obviously you need to first create a feature ' Distance to'  with the class you are interested in.  Afterwards, we going to use "find domain extreme" algorithm for finding the closest object. 

For this problem, i am going to use the same image that i used in the last blog. We are going use different concepts such as multi-threshold segmentation, variables,  and PPO. For every blob, we will find a nearest blob in the image.Subsequently, we will export results as image so that you can see what is happening.

find domain extrema
Steps:

  • multi-threshold segmentation to get  all blobs and classify as class1
  • create variable to count loop ( it will be later use to export images with distinct names)
  • Use PPO to loop through class1
    •     assign current object as curr_class
    •     update variable 
    •     find domain extreme to find nearest blob that belong to class1
    •     export image


Project to find the nearest object

Color: Red- Class 1
Color: Green:cur_class
Color :Magenta: near ( blob closest to cur-class)

gif demonstrating nearest object


Wednesday, 12 November 2014

Tutorial eCognition : Finding area of classified objects within eCognition

Problem:

I have done some classification with many classes. Now i want to know the area of individual class. Please tell a way to find it and store in a text file. 

Solution

For this kind of problem, many people would resort to finding a way in ArcGIS with exported shapes files with classification. In ArcGIS, one has to do series of operations to get the desired output (add area column and populate with area for every features, sum up area for each class). These operation are not complex but if you are not good with ArcGIS, then doing all that is a tough nut to crack.

Good thing is that there is feature in eCognition ‘area of classified objects’ that can be used for this purpose. But for that you have to create that feature for every class that you have in your project. Imagine you have more than 10 classes. Creating that features for ten different classes is boring, isn’t it? At least for me, it is cumbersome. On top of that imagine you just finished with the project you are working with and you have to perform same thing again in another project with 15 classes with different names. Gosh, you have to create another 15 features again. Not fun, right? So what’s the solution? You guess it right, we can use array functionality within eCognition to loop over all classes and create that feature for every class in one go. The rule set,we will create, will work in every project regardless of numbers of classes you have in your project. Sounds interesting? Keep on reading.

eCognition introduce concept of variables and array from version 8 ( ? not sure). With array functionality, you can gather all classes present in  you project in one go and make a particular feature with each class. I have seen many eCognition projects developed by other people for last three years and not many people use these new features of eCognition. It very handy for many cases. So here is a rule that that does the following. You can even create a customized algorithm from it and apply it to every project at the end of your classification.
  1.  create an array to store all you classes
  2.  create a temp class
  3. create a feature area of classified objects ( you can specify unit you want) based on temp class
  4.  loop each class from step 1 and store class into temp class. Get information about area from feature in step three and store in an array.
  5.  write class array and area array in a CSV with export project statistics
Area of classified objects for each class in a csv file.

Rule set in action

Loop over all classes one at a time

Store areas of current class in an array
export area information store in the array in a csv file

I have written another use of array functionality to  export a particular feature fore every image bands in your project in this blog. 


Wednesday, 5 November 2014

Python tutorial: Converting a raster dataset to XYZ in Python

One of the most popular posts in my blog is about converting a raster image into XYZ text file.  Converting raster image to XYZ file may be necessary because machine learning algorithms (outside proprietary software) requires input to be a table. I had written a post about converting a raster file into XYZ using ARCGIS some three years ago, which seems to attract many visitors to my blog.

Here I will write a way to do it using Python. For reading geo-referenced raster file, we will use rasterio package which is a wrapper for gdal that provides clean and fast I/O for geospatial raster images. The package is written in cython therefore it is very fast. Reading a raster file with rasterio is a one liner code. You can download binary of rasterio here. Its a binary file so installing rasterio is just a matter of clicking binary executable file. After reading the raster file we will get bounding box of the image and compute XY of each pixel and later write in a csv file together with pixel values. The code works for any number of bands in the image.

I hope this piece of code is useful for you. In near future I will show you to do exact same thing within QGIS. Stay tuned.


Output XYZ file with band values

In addition, the above code can be combine with the code (Data exchange between MATLAB and Python: Reading and writing .mat files with Python)  so that  you can easily export multi-spectral and hyperspectral data as a mat file for MATLAB. Many people seem to be having problem with reading remote sensing image with multiband read function in MATLAB. If you are one them, use above codes and bypass multiband read function to get your remotes images straight in MATLAB.

Tuesday, 28 October 2014

eCognition tutorial: Almost connected components for clumps identification in eCognition

This post is inspired from a post by Steve Eddins, who works in Math works, a company that build MATLAB. He is a software development manager in the MATLAB and one of the co-author of a book " Digitial image processing with MATLAB".I use both MATLAB and eCognition, so I ponder if this can be done in eCognition. eCognition has basic Morphological Operators like dilation and erosion . Advance MM operators like by opening by reconstruction, connected component labeling or skeleton and many others are not available in eCognition.

The problem of almost connected components:

There is simple synthetic image containing a number of circular blobs.  How can we label and measure the three clumps instead of the smaller circles? Two circles are almost connected if circles are within 25 pixels unit.

Binary circles
Connected components labeling
Almost connected components labeling

Of course it can be solved in eCognition. But for this, you have to  be familiar with many concepts in eCognition. Concepts such as PPO, object variables, multi-level representation and temporary layers are required. My workflow for  the solution is as follows:

  • Use multi-threshold segmentation to get circles
  • Use distance map algorithm to get binary distance map
  • Use chessboard segmentation to get pixel level unclassified objects
  • Use  multi-thresholding segmentation based on the distance map to get  clump
  • Copy level above
  • Use object variable concept to assign each clump a unique ID
  • Convert to sub-objects to get original circle at upper level

I will post rule-set after some time. I have given you enough hint how to proceed. Get your hands dirty !

Almost connected components labelling within eCognition
The concept of “almost connected components” can be applicable in remote sensing for clustering buildings detected in remote sensing images for analyzing of micro-climate of urban areas. There can be various other applications. Can you think of any ?

Monday, 27 October 2014

eCognition tutorial: Image object fusion in eCognition: Example of water bodies classification

eCognition is a powerful software for analysis for remote sensing images. Many people have this false impression that eCognition is all about segmentation. That’s not true. Segmentation is a just a part of it big chain that may involve segmentation, temporary classification, fusion, exploitation of contextual information etc. Many people believe that segmentation should be perfect at first time which is a fallacy. You can modify your segments as more information becomes available during the analysis. Typically in my any project I use segmentation at least 10 times. Yes at least 10 times .One of  the underutilized feature of eCognition is “Image Object Fusion” algorithm.The algorithm is an essential part of "iterative segmentation and classification" approach of GEOBIA. In this blog, I will show an example of image object fusion for classification of water bodies in a small remote sensing image. The step involves:
  1. Segmentation
  2. Initial classification of water
    1. Establishment of customized feature RationNIR 
      • RationNIR = (meanNIR/(meanR+meanB+meanG+meanNIR))*100
    2. Classify objects that satisfy property
      • RationNIR  less than 15 
      • Area greater than 10 pixels ( to avoid small shadows)
  3. Image object fusion using PPO ( Parent Process Object) to get whole water body
    1. Starting from water classified in step 2, check neighboring water objects and if difference between water object and neighboring objects is less than 5 in RatioNIR. Run in a loop.
    2. Process all water objects. 
Original image
Initial segmentation
Initial water classification
Final segmentation after image object fusion
Final water classification using image object fusion
Very powerful image object fusion algorithm

Step 3 is a one liner algorithm using "Image object Fusion".  Notice various parameters like Class filter, candidates classes, fitting function threshold, use absolute fitting value and weighted sum.To help you understand, I have made a "gif" to elaborate what is going on.

Blue: Initial water
Yellow: Active object
Red: Fused water after image object fusion



Here is the whole Rule set for water classification with image object fusion.

Friday, 7 February 2014

Example of extracting grid nodes in MATLAB

For last few posts, I am writing about morphological image analysis with python with open source package scikit-image. In this post I want to show benefit of morphological image analysis with some real application. I been working with road extraction from satellite images for past few weeks and  one of the common task is to extract the grid corner nodes of the road network. So let’s assume the road grid has been already extracted and we want to extract intersection of grid network. Here I am showing three approaches for that based on synthetic data. Obviously with real data, the task is more complex and might require many slight modifications.

Here is the grid. You can generate such grid with checkerboard function in MATLAB with minor processing (edge detection and hole filling). Intersections of of grid lines have a property that it has 4 white pixels in 4-N neighborhood.  We can exploit this property to extract it. I am going to show to ways by which it can be done and in addition, another method based on extraction of corners will also be shown.

%make a checkerboard
I=checkerboard(20, 10,10)>0.5;
imshow(I)
% detect edge with canny
BW = edge(I,'canny');
figure, imshow(BW)
%make square structural element to fill holes with closing
SE = strel('square', 3);
BW1 = imclose(BW, SE);
figure, imshow(BW1)
% make a skeleton of the to lines 1 pixel thick
BW2 = bwmorph(BW1,'skel',Inf);
figure, imshow(BW2)


1

 26


In mathematics and, in particular, functional analysis, convolution is a mathematical operation on two functions f and g, producing a third function that is typically viewed as a modified version of one of the original function .In simple word, convolution of image is a process where you scan your image from left to right and top to bottom within local neighborhood defined by user, and do mathematical calculation within the local neighborhood. Here we are going to use it for simply sum the number of pixels which are which are white within 3x3 neighborhood. As shown by above figure, if the point is an intersection grid point then at that point sum should be 5. Any points with less than 5 white pixels are not intersection grid points.

% with block processing and inline function with convolution
tic;
kernel = [0 1 0; ... %# Convolution kernel
1 1 1; ...
0 1 0];
sumX = conv2(double(BW2),kernel,'same');
result=sumX;
% only consider pixels which are in grid image
result (BW2==0)= 0
result(sumX<5)=0;
[r,c] = find(result>0);
toc;
figure,imshow((BW2)), hold on, plot(r,c,'r*');title ('with Convolution')
hold off
3
Way Two: With Mathematical Morphology (MM)
If you are working in the field of Geo-information and you still saying what the heck this ‘MM’ then seriously you should be acquainted with this field of image processing. You can use MM for classification, edge detection, filtering, segmentation, building detection and many other things. Grab a cup of coffee and Google MM in remote sensing; there is tons of stuff to read.  There are courses in MM taught in different universities which are one semester long, so you got an idea how vast the subject is. If you want to seriously know MM, then there is no better book than ‘ Morphological Image Analysis’ by P Soillie. Here we use ‘Erosion’ method of MM.  Time: 0.011184 seconds.

% method two with morphological analysis
tic;
SE = strel ('disk', 1);
result2 = imerode(BW, SE);
[r,c] = find(result>0);
toc;
figure,imshow((BW2)), hold on, plot(r,c,'b*'), title ('With MM Erosion')
hold off
4
Way Two:  Mathematical Morphology (MM)
There are many corner detection techniques such as Moravec, Harris, SIFT etc. Here, I am going to use Harris corner detector for intersection grid points. Time: 0.094848 seconds.

% detect corner by harris corner detector
tic;
C = corner(BW2,500 );
timet2 = toc;
toc;
figure,imshow(BW2);
hold on
plot(C(:,1), C(:,2), 'g*'); title('with harris corner detector')
hold off
5
Way Two:  With Harris Corner Detector
So you have seen there are number of ways to solve a particular problem. All above methods detected intersection points successfully. Among three methods, convolution required longest time for processing and MM required shortest time. The gain in computation cost between MM and Corner detection is 9 folds.  I did all this in MATLAB but it can easily coded in python with scikit-image.

Monday, 18 November 2013

Image processing in python with scikit-learn and pymorph

I have been using MATLAB quite a bit for some time. But being a proprietary language it has some drawbacks. One has to fork out lot of money to buy it and renew it yearly license. On top of that one hast to busy many toolboxes.Check this blog for finding various advantages of switching over from MATLAB to Python. You will find very helpful insights into key difference  between MATLAB and python numpy here.

Recently for couples of weeks now, I have to experimenting with python which is an open source software for some of my tasks that I do mainly related to classification, segmentation of remote sensing image etc during my weekend free time. So after sacrificing two of my weekends time for scouring python and different libraries which otherwise I would spend for my leisure activity, I am now pretty comfortable with python language as well for image processing tasks. Sweet :-) .

In this post I will show some morphological image analysis with python. I can see this will be series of post related to python. I plan to cover classification, clustering, edge detection in my subsequent posts.

For installing python, I have used pythonXY which comes with many libraries that might be needed for the scientific data analysis such as scipy, numpy, matplotlib and many others.PythoXY also comes with an Interactive development environment (IDE) for python which is called Spyder. With Spyder, you don’t have run python from command prompt. Other benefit includes integrated help, variable explorer, history console etc. In addition its all FREE.

Here is the interface how it looks.
Spyder iterface

I will be using library called “scikit-image” which comes in a bundle with pythonXY and in addition I would be using a library called “PyMorph” for morphological image analysis. Scikit-image also has a morphologicl module but there you would only find basic morphological operators like opening, closing, erosion and dialation. In PyMorph there are many advance morphological operator such opening by reconstruction and closing by reconstruction, Auto sequential filter(ASF), ASF by reconstruction etc. So after installing pythonXY, download the PYMoprh library and install it.

So, in this post I would show you how to use scikit-image for creating some basic shapes, calculating i area for each shapes or blobs and writing areas in the figure, finding bounding box of each blobs and plotting it in the figures. After finding areas of each figures, I will show how to remove certain blob  which has less area t than specified area.

% -*- coding: utf-8 -*-
Created on Mon Nov 18 15:42:21 2013
@author: shailesh
"""
%%%
% call necessary libraries
import math
import matplotlib.pyplot as plt
import numpy as np
from skimage.draw import ellipse
from skimage.draw import polygon
from skimage.draw import circle
from skimage.morphology import label
from skimage.measure import regionprops
from skimage.transform import rotate
import matplotlib.patches as mpatches
import pymorph as MM

%%%
% draw some arbitary shapes
image = np.zeros((1000, 1000),dtype=uint8)
% create an ellipse with centre (350,350) with minor axis = 100 and major = 220
rr, cc = ellipse(300, 350, 100, 220)
image[rr,cc] = 1

% create a polygon
x = np.array([1, 7, 4, 1])
x = np.array([1, 70, 40, 1])
y = np.array([1, 20, 80, 1])
rr, cc = polygon(y, x)
image[rr, cc] = 1

% create a polygon
rr, cc = circle(200, 200, 50)
image[rr, cc] = 1
fig, ax = plt.subplots(ncols=1, nrows=1, figsize=(6, 6))
ax.imshow(image)

%lable connected regions
label_img = label(image)
%find properties of connected regions
regions = regionprops(label_img,['Centroid', 'BoundingBox', 'Area'])

%loop through connected regions and find
% Centroid', 'BoundingBox', 'Area' for each blobs
for props in regions:
y0, x0 = props['Centroid']
minr, minc, maxr, maxc = props['BoundingBox']
area =props['Area']
rect = mpatches.Rectangle((minc, minr), maxc - minc, maxr - minr,
fill=False, edgecolor='red', linewidth=2)
ax.add_patch(rect)
ax.text(x0, y0, str(area), fontsize=10, color="blue")
%show figure
plt.gray()
plt.axis((0, 1000, 1000, 0))
plt.show()

%%% remove small blobs
image = image.astype(np.uint8)
%remove blobs with area less than 10000 '
% MM.areaopen is coming from PyMorph
b=MM.areaopen(image,10000)
fig, ax = plt.subplots(ncols=1, nrows=1, figsize=(6, 6))
ax.imshow(b)
ax.title('blobs with areas greater than 10000 ')
plt.gray()
plt.axis((0, 1000, 1000, 0))
plt.show()
% to close the pixels
plt.close("all")

Here are the figures produced by above code.

Initial figure after some blobs were made. The area of each blobs and bounding box were plotted as well
figure_1

Figure where blobs with area less than 10000 pixels are remove with pymorph.
figure_2

So this is a hypothetical simple illustration. But the real application of this kind of analysis could be many in remote sensing field. For example, after detection buildings you can easily remove buildings with size less than minimum mapping unit (say 100 m2). Similarly it can be used for removing some artifacts which are small in area that you might not be interested in. And the list can go on.

In next few blogs, i will be writing more about other morphological image processing operators as well as classfiers such Support Vector Machine (SVM), Random Forest (RF) etc.

Monday, 10 June 2013

Maximum reflectance in a spectra of Multispectral or Hyperspectral image in MATLAB

This is a very short post, something related to processing of remote sensing images. The image can be either multi-spectral or hyper-spectral images. One of my colleagues asked me this simple help in MATLAB  He was trying that in ENVI with band-math and IDL, with no avail.

I wanted to produce 2 images from my hyper-spectral image: 1) showing the maximum reflectance value across the bands and 2) the band number where the maximum came from.

First step is to read the image in MATLAB, of course :). That's no-brainer. For that you have to use multibandread function. I have a written post about it earlier. Search in the site. Now you have your image read in the MATLAB.

Here is the code, that i have written to perform what my friend asked.Hope this is useful for someone else.

%find the number of rows and columns
[rows,cols,bands]= size(imt1);
%reshape it 
imt1= reshape(imt1,[],4)';
%find max_value and band which has maximum value
[max_value,idx]= max(imt1,[],1);
%reconstruct the image
band_image= reshape(idx,rows,cols);
max_image= reshape(max_value,rows,cols);

Friday, 7 June 2013

eConition Tutorial: Customized algorithm for performing majority vote in eCognition

Today, I present you a customized rule set which lets you to assign super-object by evaluating all of its sub-objects based on which classification makes up the largest proportion of the area. This is one of the wishlist in eCognition Ideas and was also frequently asked in the ecognition community.

A majority statistic customized algorithm lets you assign  super- object to the class with the majority value of the pixels within each object. This would be useful for converting existing pixel-based classifications into an object-based format where additional object-based edits can be made.

The only parameter in this customized algorithm is a level variable. This is a level which consists of super-object. The customized algorithm will look one level below the level variable and perform all the necessary calculation.

Sub-objects level with classification


Super-objects level with no classification

Super-object level with classification from customized algorithm
The customized algorithm does assignment of super-object with following steps.

1) store all the classes in an array ( array_class)
3) loop through your objects in super objects
2) loop through your class arrays
4) store rel. area of (1) in an array for the class ( array_occur)
5) find maximum in the array_occur
6) assign super-object as maximum occurrence class in array  array_occur.

All of these steps are performed behind the scene in customized algorithm, so the user does not need to worry about how to perform these steps.

I have uploaded a zip file where you will find a project which shows the usage of the customized algorithm. There is a  customized rule-set as well. Load the customized rule set in your project and after that you will find a algorithm MajorityVote in the available algorithm list.

The algorithm was developed with eCogntion 8.8 and  will not work below that version.