package yi.or;

import java.awt.*;
import yi.iu.processor.*;
import yi.iu.processor.morph.*;
import yi.util.io.*;
import yi.iu.ui.SegFrame;

public class Region
{
    boolean VERBOSE=false;

    int[][] map= null;
    int image_width;
    int image_height;

    static int OUT=0;
    static int IN=1;

    public Region(int w, int h) {
        map= new int[h][w];
        image_width=w;
        image_height=h;
    }

    public Region(Dimension d) {
        this(d.width, d.height);
    }

    public Region(int[][] _map)
    {
        map=new int[_map.length][];
        for (int i=0;i<_map.length;i++) map[i]=(int[])_map[i].clone();
        image_height= _map.length;
        image_width= _map[0].length;
    }

    public Region(Region r)
    {
        this(r.map);
    }

    public Region(Region r, int rgn_num)
    {
        int[][] _map=r.map;
        image_height= _map.length;
        image_width= _map[0].length;
        map= new int[image_height][image_width];
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++)
                if (_map[y][x]==rgn_num) map[y][x]=IN;
    }

    public Region setAll() {
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++)
                map[y][x]=IN;
        return this;
    }

    public void setPixel(int x, int y) {
        if (map==null) System.err.println("use region before init");
        else if (x<0||x>=image_width||y<0||y>image_height)
            System.err.println("invalid coordinator for region");
        else map[y][x]= IN;
    }

    public Region setRegion(Rectangle r) {
        int left=r.x, right=left+r.width, top=r.y, bottom=top+r.height;
        for (int y=top;y<bottom;y++)
            for (int x=left;x<right;x++)
                map[y][x]=IN;
        return this;
    }

    public void unsetPixel(int x, int y) {
        if (map==null) System.err.println("use region before init");
        else if (x<0||x>=image_width||y<0||y>image_height)
            System.err.println("invalid coordinator for region");
        else map[y][x]= OUT;
    }

    public void setPixel(int index) {
        int y= index/image_width;
        int x= index- y*image_width;
        setPixel(x,y);
    }

    public void unsetPixel(int index) {
        int y= index/image_width;
        int x= index- y*image_width;
        unsetPixel(x,y);
    }

    public Region unsetAll() {
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++)
                map[y][x]=OUT;
        return this;
    }


/*******************************************************************************
morphology
*******************************************************************************/
    public Region closing(int struct, int size) {
        (new Closing(struct, size)).processOn(map);
        return this;
    }

    public Region opening(int struct, int size) {
        (new Opening(struct, size)).processOn(map);
        return this;
    }

    public Region erosion(int struct, int size) {
        (new Erosion(struct, size)).processOn(map);
        return this;
    }

    public Region dilation(int struct, int size) {
        (new Dilation(struct, size)).processOn(map);
        return this;
    }

/*******************************************************************************
connected component
*******************************************************************************/

    //relabel the region, (connected component) return the max label, the range is [0..maxLabel]
    public int connectComponent() {
        Dimension d= new Dimension(image_width, image_height);
        int[] map_1= yi.util.ImageUtil.twoD2oneD(map);
        ConnectComponent cc=new ConnectComponent();
        map_1= cc.compactLabeling(map_1, d);
        map= yi.util.ImageUtil.oneD2twoD(map_1, d);
        map_1=null; System.gc();
        return cc.getMaxLabel();
    }


/*******************************************************************************
return the result
*******************************************************************************/
    public int[] getArray()
    {
        int[] array= new int[image_width*image_height];
        int index=0;
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++) {
                array[index++]= map[y][x];
            }
        return array;
    }

    public int[] getMasked(int[] source) {
        int[] masked= (int[])source.clone();

        int index=0;
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++) {
                if (map[y][x]!=IN) masked[index]=0xff<<24;
                index++;
            }
        return masked;
    }

/*******************************************************************************
functions: filters
*******************************************************************************/
    //how to filter
    static int MAKEFALSEOUT=0;
    static int MAKETRUEIN=1;
    //give the information and the possible measure result of regions wanted
    //variations and means are given by info, in order
    public Region rangeFilter(int[][] info, int[][] range)
    {
        if (range==null) return this;
        int size=info.length;
        if (size!=image_width*image_height) {
            System.err.println("information size is different from image size!");
            return this;
        }
        int info_num= info[0].length;

        int index=0;
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++) {
                boolean ok= true;
                for (int j=0;j<info_num;j++) {
                    if (range[j]==null) continue;
                    if (info[index]==null) {
                        System.err.println("no available information on #"+index);
                        continue;
                    }
                    if ((info[index][j]<range[j][0])||(info[index][j]>=range[j][1])) {
                        ok=false;
                        break;
                    }
                }
                if (ok) map[y][x]=IN;
                index++;
            }
        return this;
    }


    //give the information and the possible measure result of regions wanted
    //variations and means are given by info, in order
    public Region distanceFilter(int[][] info, int[] expected, int dist)
    {
        if (expected==null) return this;
        int size=info.length;
        if (size!=image_width*image_height) {
            System.err.println("information size is different from image size!");
            return this;
        }
        int info_num= info[0].length;

        int index=0;
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++) {
                boolean ok= true;
                int d=0;
                for (int j=0;j<info_num;j++) {
                    int _d=info[index][j]-expected[j];
                    d+= _d*_d;
                }
                if (d<=dist) map[y][x]=IN;
                index++;
            }
        return this;
    }

    public Region rangeFilterOnRegion(double[][] info, double[][] range)
    {
        int[] region_yes= new int[info.length];
        for (int i=1;i<region_yes.length;i++) {
            boolean ok= true;
            for (int j=0;j<info[i].length;j++) {
                if (range[j]==null) continue;
                if ((info[i][j]<range[j][0])||(info[i][j]>=range[j][1])) {
                    ok=false;
                    break;
                }
            }
            if (ok) region_yes[i]=IN;
            else region_yes[i]=OUT;
        }
        region_yes[0]=OUT;

        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++) {
                if (map[y][x]>=region_yes.length) map[y][x]=OUT;
                else map[y][x]=region_yes[map[y][x]];
            }
        return this;
    }

    //give the information and the possible measure result of regions wanted
    //variations and means are given by info, in order
    public Region mvFilter(int[][] info, int[][] filter)
    {
        //change the mean-variation information to range information
        int[][] ranges= mv2range(filter,2);
        rangeFilter(info,ranges);
        return this;
    }

    //change the variation-mean requirement to minimum-maximum requirement
    int[][] mv2range(int[][] mv, double ex_ratio)
    {
        int len= mv.length;
        int[][] range= new int[len][2];
        for (int i=0;i<len;i++) {
            range[i][0]=(int)Math.round(mv[i][0]-mv[i][1]*ex_ratio);
            range[i][1]=(int)Math.round(mv[i][0]+mv[i][1]*ex_ratio);
        }
        return range;
    }


/*******************************************************************************
region inter- functions (union, intersection, difference)
*******************************************************************************/
    Region difference(Region r)
    {
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++) {
                if (map[y][x]!=IN) continue;
                if (r.map[y][x]==IN) map[y][x]=OUT;
            }
        return this;
    }

    Region intersection(Region r)
    {
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++) {
                if (map[y][x]==IN&&r.map[y][x]!=IN) map[y][x]=OUT;
            }
        return this;
    }

    Region intersection(Rectangle r)
    {
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++) {
                if (!r.contains(x,y)) map[y][x]=OUT;
            }
        return this;
    }

    Region union(Region r)
    {
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++) {
                if (map[y][x]==OUT&&r.map[y][x]!=OUT) map[y][x]=r.map[y][x];
            }
        return this;
    }

/*******************************************************************************
region properties
*******************************************************************************/
    public int getSize()
    {
        int size=0;
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++) {
                if (map[y][x]==IN) size++;
            }
        return size;
    }

    public Rectangle getExtension()
    {
        int[] hor= new int[image_width];
        int[] ver= new int[image_height];
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++) {
                if (map[y][x]==IN) {
                    hor[x]++;ver[y]++;
                }
            }
        //threshold of number of pixel to be the boundary, to eliminate the noise point
        int size=getSize()/100;
        int left,right,top,bottom;
        int sum=0;
        for (left=0;left<image_width;left++) {
            sum+=hor[left];
            if (sum>size) break;
        }
        sum=0;
        for (right=image_width-1;right>=0;right--) {
            sum+=hor[right];
            if (sum>size) break;
        }
        sum=0;
        for (top=0;top<image_height;top++) {
            sum+=ver[top];
            if (sum>size) break;
        }
        sum=0;
        for (bottom=image_height-1;bottom>=0;bottom--) {
            sum+=ver[bottom];
            if (sum>size) break;
        }
        return new Rectangle(left, top, right-left+1, bottom-top+1);
    }

    static int MAINPART=4;  //when return the range of the info, MAINPART%-(100-MAINPART) are considered
                            //eg.  MAINPART=4, 4%-96% is considered as the main part
    public int[][] getInfoStat(int[][] info)
    {
        int info_num= info[0].length;
        int[][] rst= new int[info_num][8];
        for (int i=0;i<info_num;i++)
        {
            int index=0;
            long sum=0;
            int size=0;
            int max=info[i][0], min=info[i][0];
            //the size, min, max, mean
            for (int y=0;y<image_height;y++) {
                for (int x=0;x<image_width;x++) {
                    if (map[y][x]!=OUT) {
                        sum+= info[index][i];
                        size++;
                        if (info[index][i]>max) max=info[index][i];
                        if (info[index][i]<min) min=info[index][i];
                    }
                    index++;
                }
            }
            if (size==0) continue;
            int mean=(int)(sum/size);

            //the variation, histogram
            long var=0;
            index=0;
            int[] hist= new int[max+1];
            for (int y=0;y<image_height;y++) {
                for (int x=0;x<image_width;x++) {
                    if (map[y][x]!=OUT) {
                        long diff= info[index][i]-mean;
                        var+= diff*diff;
                        hist[info[index][i]]++;
                    }
                    index++;
                }
            }
            var=(long)Math.sqrt(var/size);
            //OutputResult.outputResult("stat_hist", hist, 1, hist.length);

            //the range 2%-98%
            int threshold= (int)(size*MAINPART/100);
            int[] range= {0,0};
            if (threshold>0) range=getRange(hist,threshold);

            //record the mean and variation
            rst[i][0]=mean;
            rst[i][1]=(int)var;
            rst[i][2]=range[0];
            rst[i][3]=range[1]+1;
            rst[i][4]=min;
            rst[i][5]=max+1;
            rst[i][6]=(mean-var)<min?min:range[0];
            rst[i][7]=(mean+var)>max?max+1:range[1]+1;
        }


        int index=0;
        long sum=0;
        int size=0;
        for (int y=0;y<image_height;y++) {
            for (int x=0;x<image_width;x++) {
                int d=0;
                if (map[y][x]!=OUT) {
                    for (int i=0;i<info_num;i++) {
                        int _d= (info[index][i]-rst[i][0]);
                        d+=_d*_d;
                    }
                    sum+=d;
                    size++;
                }
                index++;
            }
        }
        System.err.println("average diff= " +sum/size);
        index=0;
        int[] stat= new int[256];
        int quan= (int)(sum/size/128);
        for (int y=0;y<image_height;y++) {
            for (int x=0;x<image_width;x++) {
                int d=0;
                if (map[y][x]!=OUT) {
                    for (int i=0;i<info_num;i++) {
                        int _d= (info[index][i]-rst[i][0]);
                        d+=_d*_d;
                    }
                    d/=quan;
                    if (d>255) d=255;
                    stat[d]++;
                }
                index++;
            }
        }


        OutputResult.outputResult("d_stat", stat, 1, 256);
        System.gc();
        return rst;
    }

    public int[] getInfoAvg(int[][] info)
    {
        int info_num= info[0].length;
        int[] rst= new int[info_num];
        for (int i=0;i<info_num;i++)
        {
            int index=0;
            long sum=0;
            int size=0;
            int max=info[i][0], min=info[i][0];
            //the size, min, max, mean
            for (int y=0;y<image_height;y++) {
                for (int x=0;x<image_width;x++) {
                    if (map[y][x]!=OUT) {
                        sum+= info[index][i];
                        size++;
                        if (info[index][i]>max) max=info[index][i];
                        if (info[index][i]<min) min=info[index][i];
                    }
                    index++;
                }
            }
            if (size==0) continue;
            int mean=(int)(sum/size);

            rst[i]=mean;
        }

        System.gc();
        return rst;
    }

    //return percentage of total image area covered by the blob
    int area(int label)
    {
        int a=0;
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++)
                if (map[y][x]==label) a++;
        return a;
    }

    //return mean x and y coordinates, expressed as fractions of the width and height of the image, respectively
    Point position(int label)
    {
        int a= area(label);
        long r=0,c=0;
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++)
                if (map[y][x]==label) {
                    r+=x;c+=y;
                }

        Point p= new Point((int)((r+a/2)/a),(int)((c+a/2)/a));
        return p;
    }

    Rectangle extent(int label, int extent)
    {
        int xs[]= new int[image_width], ys[]= new int[image_height];
        long a=0;
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++)
                if (map[y][x]==label) {
                    xs[x]++;ys[y]++;
                    a++;
                }
        int x1=0,x2=image_width-1,y1=0,y2=image_height-1;
        int threshold= (int)Math.round(a*extent/100);
        int sum=0;
        for (int x=0;x<image_width;x++) {
            sum+=xs[x];
            if (sum>=threshold) {
                x1=x; break;
            }
        }
        sum=0;
        for (int x=image_width-1;x>=0;x--) {
            sum+=xs[x];
            if (sum>=threshold) {
                x2=x; break;
            }
        }
        sum=0;
        for (int y=0;y<image_height;y++) {
            sum+=ys[y];
            if (sum>=threshold) {
                y1=y; break;
            }
        }
        sum=0;
        for (int y=image_height-1;y>=0;y--) {
            sum+=ys[y];
            if (sum>=threshold) {
                y2=y; break;
            }
        }

        Rectangle r= new Rectangle(x1,y1,x2-x1,y2-y1);
        return r;
    }

    static int EXTENT=5;       //20%-80% is considered as the extent
    Rectangle extent(int label)
    {
        return extent(label, EXTENT);
    }

    double aspect(int label)
    {
        Rectangle r= extent(label);
        double dx= r.width;
        double dy= r.height;
        return Math.abs((dx-dy)/Math.sqrt(dx*dy));
    }

    double slant(int label)
    {
        Point p= position(label);
        long d1=0, d2=0;
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++)
                if (map[y][x]==label) {
                    int t=(x-p.x)*(y-p.y);
                    d1+=t;
                    d2+=Math.abs(t);
                }
        return Math.abs((double)d1/d2);
    }

    double density(int label)
    {
        Rectangle r= extent(label);
        int a=0;
        for (int y=r.y; y<r.y+r.height; y++)
            for (int x=r.x; x<r.x+r.width; x++)
                if (map[y][x]==label) {
                    a++;
                }
        return (double)a/(r.width*r.height);
    }


    static int CONNECTED4=4;    //the returned perimeter is 4-connected
    static int CONNECTED8=8;    //the returned perimeter is 8-connected
    int[][] perimeter(int label, int connectivity, boolean frameyes)
    {
        int[][] p= new int[image_height][image_width];
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++) {
                if (map[y][x]!=label) continue;
                boolean isP=false;
                if ((y>0)?map[y-1][x]!=label:frameyes) isP=true;
                else if ((y<image_height-1)?map[y+1][x]!=label:frameyes) isP=true;
                else if ((x>0)?map[y][x-1]!=label:frameyes) isP=true;
                else if ((x<image_width-1)?map[y][x+1]!=label:frameyes) isP=true;
                else if (connectivity==CONNECTED4) {
                    if ((y>0&&x>0)?map[y-1][x-1]!=label:frameyes) isP=true;
                    else if ((y<image_height-1&&x>0)?map[y+1][x-1]!=label:frameyes) isP=true;
                    else if ((y>0&&x<image_width-1)?map[y-1][x+1]!=label:frameyes) isP=true;
                    else if ((y<image_height-1&&x<image_width-1)?map[y+1][x+1]!=label:frameyes) isP=true;
                }
                if (isP) p[y][x]=label;
            }
        return p;
    }

    //calculate the perimeter length in the given edge bitmap
    double perimeterLength(int[][] p, int label)
    {
        double sqrt2= Math.sqrt(2);
        Region rgn= new Region(p);
        double len=0;
        while (true) {
            //find a start point
            int bx=0, by=0;
            boolean found=false;
            for (int y=0;y<image_height;y++)
                for (int x=0;x<image_width;x++) {
                    if (p[y][x]==label) {
                        by=y; y=image_height;
                        bx=x;
                        found=true;
                        break;
                    }
                }
            if (!found) break;
            while (true)
            {
                if ((by>0)?p[by-1][bx]==label:false) {
                    p[by][bx]=OUT; by--;
                    len++; continue;
                }
                else if ((by<image_height-1)?p[by+1][bx]==label:false){
                    p[by][bx]=OUT; by++;
                    len++; continue;
                }
                else if ((bx>0)?p[by][bx-1]==label:false){
                    p[by][bx]=OUT; bx--;
                    len++; continue;
                }
                else if ((bx<image_width-1)?p[by][bx+1]==label:false){
                    p[by][bx]=OUT; bx++;
                    len++; continue;
                }
                else if ((by>0&&bx>0)?p[by-1][bx-1]==label:false) {
                    p[by][bx]=OUT; by--; bx--;
                    len+=sqrt2; continue;
                }
                else if ((by<image_height-1&&bx>0)?p[by+1][bx-1]==label:false) {
                    p[by][bx]=OUT; by++;bx--;
                    len+=sqrt2; continue;
                }
                else if ((by>0&&bx<image_width-1)?p[by-1][bx+1]==label:false) {
                    p[by][bx]=OUT; by--;bx++;
                    len+=sqrt2; continue;
                }
                else if ((by<image_height-1&&bx<image_width-1)?p[by+1][bx+1]==label:false) {
                    p[by][bx]=OUT; by++;bx++;
                    len+=sqrt2; continue;
                }
                p[by][bx]=OUT;
                break;
            }
        }
        return len;
    }
    double perimeterLength(int label)
    {
        int[][] p= perimeter(label,CONNECTED8,true);
        return perimeterLength(p, label);
/*
        double sqrt2= Math.sqrt(2);
        Region rgn= new Region(p);
        int area= rgn.area(label);
        double len=0;
        while (area>0) {
            //find a start point
            int bx=0, by=0;
            for (int y=0;y<image_height;y++)
                for (int x=0;x<image_width;x++) {
                    if (p[y][x]==label) {
                        by=y; y=image_height;
                        bx=x;
                        break;
                    }
                }
            while (true)
            {
                if ((by>0)?p[by-1][bx]==label:false) {
                    p[by][bx]=OUT; by--;
                    area--;
                    len++; continue;
                }
                else if ((by<image_height-1)?p[by+1][bx]==label:false){
                    p[by][bx]=OUT; by++;
                    area--;
                    len++; continue;
                }
                else if ((bx>0)?p[by][bx-1]==label:false){
                    p[by][bx]=OUT; bx--;
                    area--;
                    len++; continue;
                }
                else if ((bx<image_width-1)?p[by][bx+1]==label:false){
                    p[by][bx]=OUT; bx++;
                    area--;
                    len++; continue;
                }
                else if ((by>0&&bx>0)?p[by-1][bx-1]==label:false) {
                    p[by][bx]=OUT; by--; bx--;
                    area--;
                    len+=sqrt2; continue;
                }
                else if ((by<image_height-1&&bx>0)?p[by+1][bx-1]==label:false) {
                    p[by][bx]=OUT; by++;bx--;
                    area--;
                    len+=sqrt2; continue;
                }
                else if ((by>0&&bx<image_width-1)?p[by-1][bx+1]==label:false) {
                    p[by][bx]=OUT; by--;bx++;
                    area--;
                    len+=sqrt2; continue;
                }
                else if ((by<image_height-1&&bx<image_width-1)?p[by+1][bx+1]==label:false) {
                    p[by][bx]=OUT; by++;bx++;
                    area--;
                    len+=sqrt2; continue;
                }
                p[by][bx]=OUT;
                area--;
                break;
            }
        }
        return len;
*/
    }

    double compactness(int label)
    {
        double pl= perimeterLength(label);
        double a= area(label);
        return ((pl/4)*(pl/4)-a)/a;
    }

    int getMaxLabel()
    {
        int max=0;
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++)
                if (map[y][x]>max) max=map[y][x];
        return max;
    }

    public Region removeSmallRegion(int size)
    {
        int max= connectComponent();
        int[] a= areas(max);
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++)
                if (a[map[y][x]]<size) map[y][x]=OUT;
        return this;
    }

    int relabel()
    {
        //label first
        int maxLabel= getMaxLabel();
        int[] stat= new int[maxLabel+1];
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++)
                {
                    if (map[y][x]>maxLabel)
                        System.err.println("bigger label than maxLabel found!");
                    stat[map[y][x]]++;
                }

        stat[0]=0;
        int j=1;
        for (int i=1;i<stat.length;i++) {
            if (stat[i]!=0) stat[i]=j++;
        }
        if (VERBOSE) System.err.println("From "+maxLabel+" to "+(j-1)+" regions");
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++)
                map[y][x]= stat[map[y][x]];
        return (j-1);

    }

/******************************************************************************
get the static information of all the connected component
******************************************************************************/
    int[] areas(int maxLabel)
    {
        int[] a=new int[maxLabel+1];
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++)
                a[map[y][x]]++;
        return a;
    }

    Point[] positions(int maxLabel)
    {
        int[] a=areas(maxLabel);
        long[] r=new long[maxLabel+1],c=new long[maxLabel+1];
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++) {
                int label=map[y][x];
                r[label]+=x;
                c[label]+=y;
            }
        Point[] p= new Point[maxLabel+1];
        for (int i=0;i<=maxLabel;i++) {
            p[i]= new Point((int)((r[i]+a[i]/2)/a[i]),(int)((c[i]+a[i]/2)/a[i]));
        }

        return p;
    }

    //return the min and max of the range which are the boundary of main part of whole data
    //for example, can return the range from 100-900, if limit=100, and total =1000;
    int[] getRange(int[] hist, int number)
    {
        int len= hist.length;
        int sum=0;
        int l,r;
        for (l=0;l<len;l++) {
            sum+=hist[l];
            if (sum>=number) break;
        }
        sum=0;
        for (r=len-1;r>=0;r--) {
            sum+=hist[r];
            if (sum>=number) break;
        }
        int[] range={l,r};
        return range;
    }

    Rectangle[] extents(int maxLabel, int extent)
    {
        int[][] xs= new int[maxLabel+1][image_width], ys= new int[maxLabel+1][image_height];
        int[] a= new int[maxLabel+1];
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++) {
                int label=map[y][x];
                xs[label][x]++;ys[label][y]++;
                a[label]++;
            }

        Rectangle[] r= new Rectangle[maxLabel+1];
        for (int i=0;i<maxLabel+1;i++) {
            int threshold= (int)Math.round(a[i]*extent/100);
            int[] rx=getRange(xs[i],threshold);
            int[] ry=getRange(ys[i],threshold);
            r[i]= new Rectangle (rx[0],ry[0],rx[1]-rx[0],ry[1]-ry[0]);
        }
        return r;
    }

    Rectangle[] extents(int maxLabel)
    {
        return extents(maxLabel, EXTENT);
    }

    double[] aspects(int maxLabel)
    {
        Rectangle[] r= extents(maxLabel);
        double[] as= new double[maxLabel+1];
        for (int i=0;i<as.length;i++) {
            double dx= r[i].width;
            double dy= r[i].height;
            as[i]=Math.abs((dx-dy)/Math.sqrt(dx*dy));
        }
        return as;
    }

    double[] slants(int maxLabel)
    {
        Point[] p= positions(maxLabel);
        long[] d1=new long[maxLabel+1], d2=new long[maxLabel+1];
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++) {
                int label=map[y][x];
                int t=(x-p[label].x)*(y-p[label].y);
                d1[label]+=t;
                d2[label]+=Math.abs(t);
            }
        double[] sl= new double[maxLabel+1];
        for (int i=0;i<=maxLabel;i++) sl[i]=Math.abs((double)d1[i]/d2[i]);
        return sl;
    }

    double[] densities(int maxLabel)
    {
        Rectangle[] r= extents(maxLabel);
        int[] a=new int[maxLabel+1];
        for (int i=0;i<=maxLabel;i++)
            for (int y=r[i].y; y<r[i].y+r[i].height; y++)
                for (int x=r[i].x; x<r[i].x+r[i].width; x++)
                    if (map[y][x]==i) a[i]++;
        double[] d= new double[maxLabel+1];
        for (int i=0;i<=maxLabel;i++)  d[i]=(double)a[i]/(r[i].width*r[i].height);
        return d;
    }


    //connectivity: the returned perimeter will be 4-connected (CONNECTED4) or 8-connected (CONNECTED8)
    //frameyes: whether the bounday of the image is considered as the perimeter
    int[][] perimeters(int maxLabel, int connectivity, boolean frameyes)
    {
        int[][] p= new int[image_height][image_width];
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++) {
                int label=map[y][x];
                boolean isP=false;
                if ((y>0)?map[y-1][x]!=label:frameyes) isP=true;
                else if ((y<image_height-1)?map[y+1][x]!=label:frameyes) isP=true;
                else if ((x>0)?map[y][x-1]!=label:frameyes) isP=true;
                else if ((x<image_width-1)?map[y][x+1]!=label:frameyes) isP=true;
                else if (connectivity==CONNECTED4) {
                    if ((y>0&&x>0)?map[y-1][x-1]!=label:frameyes) isP=true;
                    else if ((y<image_height-1&&x>0)?map[y+1][x-1]!=label:frameyes) isP=true;
                    else if ((y>0&&x<image_width-1)?map[y-1][x+1]!=label:frameyes) isP=true;
                    else if ((y<image_height-1&&x<image_width-1)?map[y+1][x+1]!=label:frameyes) isP=true;
                }
                if (isP) p[y][x]=label;
            }
        return p;
    }

    int[][] perimeters(int connectivity, boolean frameyes)
    {
        int maxLabel= getMaxLabel();
        return perimeters(maxLabel, connectivity, frameyes);
    }

    double[] perimeterLengths(int maxLabel)
    {
        int[][] ps= perimeters(maxLabel, CONNECTED8, true);
        double[] pls= new double[maxLabel+1];
        for (int i=1;i<=maxLabel;i++) pls[i]=perimeterLength(ps, i);
        pls[0]=0;
        return pls;
    }

    double[] compactnesses(int maxLabel)
    {
        double[] pl= perimeterLengths(maxLabel);
        int[] a= areas(maxLabel);
        double[] com= new double[maxLabel+1];
        for (int i=1;i<=maxLabel;i++)  com[i]=((pl[i]/4)*(pl[i]/4)-a[i])/a[i];
        com[0]=0;
        return com;
    }

    public double[][] getStatistics(int maxLabel)
    {
        double[][] stat= new double[maxLabel+1][12];

        int[] a= areas(maxLabel);
        Point[] p= positions(maxLabel);
        Rectangle[] r= extents(maxLabel);
        double[] as= aspects(maxLabel);
        double[] s= slants(maxLabel);
        double[] d= densities(maxLabel);
        double[] pe= perimeterLengths(maxLabel);
        double[] c= compactnesses(maxLabel);

        for (int i=0;i<=maxLabel;i++) {
            //area
            stat[i][0]=(double)a[i]/(image_width*image_height);
            //position
            stat[i][1]=(double)p[i].x/image_width;
            stat[i][2]=(double)p[i].y/image_height;
            //extent
            stat[i][3]=(double)r[i].x/image_width;
            stat[i][4]=(double)r[i].y/image_height;
            stat[i][5]=(double)r[i].width/image_width;
            stat[i][6]=(double)r[i].height/image_height;
            //aspect
            stat[i][7]=as[i];
            //slant
            stat[i][8]=s[i];
            //density
            stat[i][9]=d[i];
            //compactness
            stat[i][10]=c[i];
            //perimeter
            stat[i][11]=pe[i];
        }
        return stat;
    }

    double[][] getStatistics()
    {
        int maxLabel= getMaxLabel();
        return getStatistics(maxLabel);
    }

/*******************************************************************************
region display
*******************************************************************************/
    public void display()
    {
        int[] map_1= yi.util.ImageUtil.twoD2oneD(map);
        int[] img= yi.util.ImageUtil.makeNewImageFromLabel(map_1);
        addPerimeters(img);
        (new SegFrame("regions", img, image_width, image_height, map_1)).show();
        map_1=null; img=null;
        System.gc();
    }

    public void display(int[] source)
    {
        display("region",source);
    }

    public void display(String name, int[] source)
    {
        int[] img= (int[])source.clone();
        addPerimeters(img);
        int[] map_1= yi.util.ImageUtil.twoD2oneD(map);
        (new SegFrame(name, img, image_width, image_height, map_1)).show();
        map_1=null; img=null;
        System.gc();
    }
/*******************************************************************************
region other utility
*******************************************************************************/
    //in some application, the region #0 means that region belonging to no region
    //in others, region with number0 is a actual region
    //in some case, region with #0 need to be reassign a number so that it is considered as a region
    public void enable0()
    {
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++)
                map[y][x]++;
    }

    public void changeLabel(int i, int j)
    {
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++)
                if (map[y][x]==i) map[y][x]=j;
    }

    //sometimes, the labels are not from 0 to N.  there are some holes in the number of labels
    //make the labels from 0 to N continuously
    //return the max label now
    public int recalLabel()
    {
        int maxLabel=getMaxLabel();
        int[] count= areas(maxLabel);
        count[0]=0;
        int new_index=1;
        for (int i=1;i<count.length;i++) {
            if (count[i]!=0) count[i]=new_index++;
            else count[i]=-1;
        }
        for (int y=0;y<image_height;y++)
            for (int x=0;x<image_width;x++)
                map[y][x]=count[map[y][x]];
        return new_index-1;
    }

    public void addPerimeters(int[] img)
    {
        int[][] ps= this.perimeters(CONNECTED8, true);
        int[] ps_1= yi.util.ImageUtil.twoD2oneD(ps);
        int[] ps_img= yi.util.ImageUtil.makeNewImageFromLabel(ps_1);
        for (int i=0;i<img.length;i++)
            if (ps_1[i]!=0) img[i]=ps_img[i];
        ps=null;
        ps_1=ps_img=null;
        System.gc();
    }
}