1.一种区域灰度和方向成分算子的遥感影像中树冠直径提取方法,包括以下步骤:S1,输入由红、绿、蓝三个波段构成的遥感影像Image,将Image转换为256级灰度的灰度影像GrayImage,建立树木直径探测结构RStruct;输入区域变化容忍度最大阈值MaxRate、区域变化容忍度最小阈值MinRate,阈值阶梯Gap; 指定当前容忍度变量CurrentRate=MaxRate;
S101,输入遥感影像Image;
S102,将Image转换为256级灰度的灰度影像GrayImage;
S103,对GrayImage的所有像元进行统计,获得GrayImage的所有像元的均值IMean;
S104,建立树木直径探测结构RStruct,RStruct对应GrayImage上的一个位置,该结构包含以下属性;
RX,该结构在GrayImage上对应的横坐标;
RY,该结构在GrayImage上对应的纵坐标;
RCircleVector圆区域灰度和方向成分矢量,为一个8个元素的矢量,初始值为[0,0,0,
0,0,0,0,0];
RStoped标志是否需要进行继续计算,值为TRUE或FALSE,初始值为FALSE;
RScale圆形区域尺度,初始值为2;
RRemoved删除标记,值为TRUE或FALSE,初始值为FALSE;
S105,输入区域变化容忍度最大阈值MaxRate、区域变化容忍度最小阈值MinRate,阈值阶梯Gap;MaxRate的默认值为0.5,MinRate的默认值为0.1,Gap的默认值为0.01;
S106,CurrentRate=MaxRate;
S2,建立圆区域灰度和方向成分算子CircleOperator,该算子的输入为横坐标CX,纵坐标CY,尺度CScale,输出为圆区域灰度和方向成分矢量CircleVector;
S201,建立圆区域成分算子CircleOperator,该算子的输入为横坐标CX,纵坐标CY,尺度CScale,
S202, pixels=在GrayImage上取出以(CX, CY)为中心以CScale为半径的圆的范围内的所有像元;
S203,pixelnum= pixels内像元的个数;
S204,基于像元灰度值对pixels的所有像元进行K‑Means算法聚类,聚类中心点个数为
4,K‑Means初始化中心点的灰度值分别为0、64、128、182,聚类的结果将pixels中的所有像元分配到4个像元分组内,分别为group1、group2、group3、group4;
S205,计算如下灰度成分变量,a1=group1中像元个数/ pixelnum,a2=group2中像元个数/ pixelnum,a3=group3中像元个数/ pixelnum,a4=group4中像元个数/ pixelnum;
S206,初始化方向成分数组b,该数组包含4个元素,初始值为0;
S207,像元分组加和值sumvalue=0;
S208,像元计数器pixelcounter=1;
S209,tpixel=取出pixels中的第pixelcounter个元素,获取tpixel在GrayImage上的横坐标TCX和纵坐标TCY;
S210,查找tpixel所处的像元分组,如果在group1内则tvalue=1,如果在group2内则tvalue=2,如果在group3内则tvalue=3,如果在group4内则tvalue=4;
S211,建立向量v1, 该向量的起点为(CX,CY) 终点为(CX+1,CY);建立向量v2,该向量的起点为(CX,CY)终点为(TCX,TCY);计算以v1的方向为起始,以(CX,CY)为中心逆时针旋转到v2的方向对应的夹角angle;
S212,计算位置变量position=round(angle/90),其中round为进行四舍五入;
S213, 设置方向成分数组b的第position个元素的值,b[position]= b[position]+tvalue
S214,sumvalue=sumvalue+tvalue;
S215,pixelcounter=pixelcounter+1;
S216,如果pixelcounter>pixelnum则转到S217,否则转到S209;
S217,计算方向成分数组b,对应公式为b[1]=b[1]/sumvalue, b[2]=b[2]/sumvalue, b[3]=b[3]/sumvalue, b[4]=b[4]/sumvalue;
S218, 建立圆区域灰度和方向成分矢量CircleVector=[a1,a2,a3,a4,b[1],b[2],b[3],b[4]]
S219, 将CircleVector作为结果返回;
S3,对于GrayImage中的每一个像元,建立树木直径探测结构RStruct,放入RStruct列表RStructList中,RStructList的元素个数为RNum;
S301,建立RStruct列表RStructList,RStructList初始状态下包含0个元素;
S302,CTPixel=在GrayImage中取出一个像元;
S303,CTRStruct=建立一个RStruct结构;
S304,CTRStruct.RX=CTPixel在GrayImage中的横坐标,CTRStruct.RY=CTPixel在GrayImage中的纵坐标;
S305,CTRStruct.RCircleVector=CircleOperator(CTRStruct.RX,CTRStruct.RY,CTRStruct.RScale),利用圆区域灰度和方向成分算子计算CTRStruct对应位置的圆区域灰度和方向成分矢量;
S306,将CTRStruct加入到RStructList中;
S307,如果GrayImage还有项目没有处理那么转到S302,否则转到S308;
S308,RNum=RStructList的元素个数;
S4,建立非树木与重合区域删除算子DeleteOperator,该算子的输入为一个在列表中的位置TX,判断RStructList的第TX个元素是否需要删除,如果需要则将第TX个元素的RRemoved属性标记为TRUE;
S401,建立非树木与重合区域删除算子DeleteOperator,该算子的输入为一个在列表中的位置TX;
S402,DORStruct=取出RStructList中第TX个元素;
S403,DOPixels=在GrayImage上取出以(DORStruct.RX, DORStruct.RY)为中心,以DORStruct.RScale为半径的所有像元;
S404,删除比对计数器DOCounter=1,初始像元数量DOPNum=统计DOPixels中的像元个数;
S405,如果DOCounter与TX相等那么转到S410,否则转到S406;
S406,DORStruct2=取出RStructList中第DOCounter个元素;
S407,如果DORStruct2.RRemoved等于TRUE则转到S410,否则转到S408;
S408,DOPixels2=在GrayImage上取出以(DORStruct2.RX, DORStruct2.RY)为中心,以DORStruct2.RScale为半径的所有像元;
S409,基于像元位置计算DOPixels和DOPixels2的交集union,将union的所有像元从DOPixels中删除;
S410,DOCounter= DOCounter+1;
S411, 如果DOCounter<=RNum则转到S405,否则转到S412;
S412,DOResult=统计DOPixels的元素个数/DOPNum;
S413, 如果DOResult<0.75则转到S414,否则转到S415;
S414,DORStruct.RRmoved=TRUE,DORStruct.RStoped=TRUE;
S415, 第S4步骤处理结束;
S5,建立区域扩大算子EnlargeOperator,该算子的输入为一个在列表中的位置EX,判断RStructList的第EX个元素是否扩大其对应区域,如果不可以则将第EX个元素的RStoped属性标记为TRUE;
S501, 建立区域扩大算子EnlargeOperator,该算子的输入为一个在列表中的位置EX;
S502,EORstruct=取出RStructList中的第EX个元素;
S503,vector1= EORstruct.RCircleVector;
S504,vector2=CircleOperator(EORstruct.RX,EORstruct.RY, EORstruct.RScale+
1),计算EORstruct的RScale扩大1之后的圆区域灰度和方向成分矢量;
S505,EOResult=l2norm(vector1‑vector2)/l2norm(vector1),其中l2norm为l2范数;
S506,如果EOResult
S507,EORstruct.RSCale= EORstruct.RSCale+1, 转到S509;
S508,EORstruct.RStoped=TRUE;
S509, 第S5步骤处理结束;
S6,利用非树木与重合区域删除算子DeleteOperator和区域扩大算子EnlargeOperator处理RStructList,输出树冠直径检测结果;
S601,处理计数器DTCounter=1;
S602,DTRstruct=取出RStructList中的第DTCounter个元素;
S603, 如果DTRstruct.RRemoved=FALSE则执行DeleteOperator(DTCounter);
S604, 如果DTRstruct.RStoped=FALSE则执行EnlargeOperator(DTCounter);
S605, DTCounter=DTCounter+1;
S606, 如果DTCounter>RNum则转到S607,否则转到S602;
S607, CurrentRate=CurrentRate‑Gap;
S608, 如果CurrentRate
S609,如果RStructList包含的元素中存在RStoped属性为FALSE的元素,则转到S601;
否则转到S610;
S610,DTCounter=1;
S611,DTRstruct2=取出RStructList中的第DTCounter个元素;
S612, 如果DTRstruct.RRemoved=FALSE则转到S613,否则转到S614;
S613,检测到树冠,以 (DTRstruct2.RX, DTRstruct2.RY)为中心,以DTRstruct2.RScalse为半径绘制圆形,并输出该圆形的直径为2×DTRstruct2.RScale;
S614,DTCounter= DTCounter+1;
S615,如果DTCounter>RNum则转到S616,否则转到S611;
S616,第S6步骤处理结束。