Jian Wang
email: wangjian@njau.edu.cn

Hongbing Yang
email: hbyang@njau.edu.cn

Research on effect of ventilated region to grape rain shed performance based on CFD

Keywords:

Two different kinds of multi-span grape rain sheds are constructed and the wind field and rainfall field distribution are simulated based on CFD and computer simulation techniques. The airflow distribution, temperature variation, wind pressure distribution and rainfall vectors distribution are compared under two conditions with ventilated region and without ventilated area in multi-span grape rain shed. The research results show that set the ventilated area on the top of grape rain shed can improve effectively the high temperature and high humidity climate environment, and it has a benefit to grapes fertility.

Introduction

Rainproof cultivation is one important form of grape protected cultivation in southern China. Grape rainproof cultivation has two main merits. On the one hand, the grape fruit set is promoted, and the grape quality is improved. On the other hand, the soil moisture is controlled effectively, and the fungus spread and the amount of pesticide spraying are decreased significantly [START_REF] Li | Effect of rainproof cultivation on grape disease incidence and light intensity under the shelter[END_REF]. However, under the effects of the temperature rise rapidly and grape transpiration, the temperature and humidity are very high on the top of grape rain shed during May to July, and it can be damaged the grape fruit set. If the ventilated region is set on the top of shelter shed, the high temperature and humidity conditions may be improved. Nevertheless, the construction of grape rain shed is changed when ventilated area is set, and the climate environments of shelter shed how to change when ventilated region is set. These new problems need to be researched in detail. In recent years, the computer simulation techniques have a wide range of applications in agricultural facilities. Briassoulis and Mistriotis [START_REF] Briassoulis | Integrated structural design methodology for agricultural protecting structures covered with nets[END_REF] proposed the integrated structural design methodology for agricultural protecting structures covered with nets based on computer simulation method. Harmanto [START_REF] Harmanto | Microclimate and air exchange rates in greenhouses covered with different nets in the humid tropics[END_REF] simulated the airflow distribution of greenhouse under different mesh sizes of nets and predicted the air exchange rate. M. Teitela [START_REF] Teitel | Effect of wind direction on greenhouse ventilation rate, airflow patterns and temperature distributions[END_REF]researched the effect of  2 Corresponding author wind direction on greenhouse ventilation rate, airflow patterns and temperature distributions based CFD.

The aims of this research were to analysis the effect of ventilated region to climate environment, and evaluate the possibility of installing ventilated region in the grape rain shed. The airflow distribution, temperature variation, wind pressure distribution and water vectors distribution are compared under two conditions of ventilated region and non-ventilated area in grape rain shed.

2

Material and method

Grape Multi-span greenhouse

The grape multi-span greenhouse was located in the coastal areas of Jiangsu province. Wind pressure formula [START_REF] Zhou | Modern Greenhouse Engineering[END_REF] is give by

2 0 1 2 v ao p C v   (1)
Where C is pressure shape factor, ao  is outdoor air density( 3 kg/m )and 0 v is outdoor wind speed ( m/s ).

Boulard and Baille [START_REF] Boulard | Modeling of air exchange rate in greenhouse equipped with continuous roof vents[END_REF][START_REF] Boulard | Natural Ventilation Performance of Six Greenhouse and Tunnel Types[END_REF] put forward the ventilation rate formula. When only considering the role of wind pressure, the ventilation rate G is

0.5 () 2 w S G AC v  ( 2
)
Where S is vent cross-sectional area( 2 m ),A is vent distribution coefficient, w C is wind pressure effect coefficient and v is average wind speed( m/s ).

CFD numerical simulation method

Computational fluid dynamics (CFD) is a sophisticated design and analysis tool that allows us to study different design limitations based on a computational model. The CFD software used in the present work was FLUENT 6.2. Its code uses a discretization procedure of finite volumes and it incorporates the equations that govern heat flow and transfer, such as the continuity equation, momentum conservation and energy conservation. The overall methodology used to carry out a complete study with the commercial software used in this work is summarized as follows: generating the geometry of the fluid domain, meshing, defining the boundary conditions and the physical properties of the model (pre-processing), solving and finally analysis (post-processing) [START_REF] Franco | Aerodynamic analysis and CFD simulation of several cooling pads used in Mediterranean greenhouses[END_REF]. The classical mass, momentum, energy and concentration equations can be represented for a steady-state, three-dimensional flow with the following conservation equation.

  j j j j uS t x x x                    (3)
Where Φ stands for the variables of interest, i.e. the three velocity components μ j (ms -1 ), the temperature T (K), Г Φ and S Φ represent the diffusion coefficient and source term of Φ and a description of their forms can be found in Boulard [START_REF] Bartzanas | Effect of Vent Arrangement on Windward Ventilation of a Tunnel Greenhouse[END_REF][START_REF] Boularda | Computer fluid dynamics prediction of climate and fungal spore transfer in a rose greenhouse[END_REF]. The system of equations built with these variables is numerically solved with the finite volume method. Algorithms and methods for the resolution of this system of equations can be found in the fluent 6.2 user's guide. To model the turbulent constraints using the standard k-ε turbulence model , in Eq. 1 Φ stands also for the turbulent kinetic energy k (m 2 s -2 ) and the dissipation of the turbulent kinetic energy ε (m 2 s -3 ).

The airflow field and rain field of multi-span grape rain shed are researched differently. The airflow field distribution is simulated when the wind direction is 0 degree, wind velocity is 20m/s and temperature is 308K. The rain field distribution is researched when rain direction 0 degree, rain velocity is and 5m/s and temperature is 308K.

3

Results and discussion

Airflow distribution simulation of multi-span grape rain shed based on CFD

Wind pressure distribution contours of two kinds of multi-span grape rain shed are shown in figure 2. From figure 2 we can see that outer surface of two different kinds of shelter sheds are lift state under the action of wind pressure, and their external flow field distribution is almost the same. However, the internal flow field of two kinds of rain shelters is different. The flow field distribution of span 1 and 2 with ventilated region is uniform, and the flow field distribution of span 3 and 4 with ventilated region is more uneven than non-ventilated region.

Fig. 2. Pressure distribution contours of two kinds of grape rain sheds

The wind pressure distribution on the top of shelter shed is given in figure 3. The abscissa axis represents the span of grape rain shed, and the ordinate axis stands for the wind pressure distribution. It can be seen that the wind pressure distribution of Span 1 and span 2 is similar to two different kinds of grape rain shed. The wind pressure distribution of span 3 and span 4 is different to two different kinds of grape rain shed. From figure 3(a) we can know that the top of shelter shed is lift mainly. The wind lift maximal is 370 Pa and wind pressure maximal is 72Pa in the span 1and span 2, and from Span 3 to span 4, the value range of wind of upward lift is 10~189Pa. It can be seen from figure 3(b) that wind lift maximal is 354Pa and wind pressure maximal is 216Pa in the span 1and span 2, and from Span 3 to span 4, the value range of wind pressure is 0~74Pa. Simulation results show that the wind pressure value is decreased on the top of grape rain shed when the ventilated region is installed. According to formula (1)we can see that it is proportional relationship between wind speed and air pressure, So the wind speed is decreased, the wind pressure is also decreased. The airflow velocity vectors distribution of multi-span grape rain shed with ventilated region is shown in figure 4.When the ventilated region is set, the velocity of external grape rain shed is decreased, so the wind pressure value is decreased correspondingly. decreased significantly. This is consistent with formula (2), and when the internal wind speed is increased, the ventilation rate is increased accordingly, so temperature and humidity of grape rain shed with ventilated region are deceased clearly.The rainproof performance has some impact when the ventilated region is set in grape rain shed. When rainfall vertical ground and the velocity is 5m/s, the rainproof effect of two kinds of grape rain shed is simulated in figure 6. Compared figure 6(a) and figure 6(b), we can see that the rainproof performance without ventilated region is better than grape rain shed with ventilated region. However, it can be seen that rainfall distribution mainly at the near ground in grape rain shed clearly and these rainfall has little effect on the grape growth.

Summary

The two kinds of models of grape rain shed are constructed, and the airflow distribution, temperature variation, wind pressure distribution and rainfall vectors distribution are compared under two conditions of ventilated region and nonventilated area in grape rain shed based on computer simulation. The research results show that set the ventilated area on the top of grape rain shed can improve effectively the high temperature and high humidity climate environment, and it has a benefit to grapes fertility. Of course, the construction costs of grape rain shed is increased when installing ventilated region.

Fig. 1 .

 1 Fig. 1. Schematic diagram multi-span grape rain shed (mm)

Fig. 3 .Fig. 4 .Fig. 5 .Fig. 6 .

 3456 Fig. 3. Pressure distribution curve on the top of grape rain shed

Acknowledgments. This work has been funded by Doctoral Fund of Ministry of Education (No.20090097120014).


Hu Zhao
email: zhaohu@whu.edu.cn

Zhengwei Yang
email: zhengwei_yang@nass.usda.gov

Liping Di

Zhiyuan Pei
email: peizhiyuan@tom.com

Evaluation of temporal resolution effect in remote sensing based crop phenology detection studies

Keywords: Crop phenology, temporal resolution evaluation, least square, double logistic function fitting

Remote sensing based phenology detection method has been employed to study agriculture, forestry and other vegetations for its potential to reflect the variations in climate change. These studies usually utilized time series Normalized Difference Vegetation Index (NDVI) generated from various sensors through a Maximum Value Compositing (MVC) process, which minimized the contamination from cloud and simultaneously introduce degradation of temporal accuracy. In this study, we assess the impact of temporal resolution on crop phenology derivation researches by comparing three different Moderate Resolution Imaging Spectroradiometer (MODIS) datasets: daily surface reflectance, 8 day composited surface reflectance and 16 day composited NDVI. The surface reflectance data were first filtered by employing auxiliary data which contained quality and viewing geometry information, and then used to calculate NDVI with specific date. A least square method was taken to fit the survival data points to double logistic function. And finally, seven time-related metrics were obtained and matched with field observation crop phenology stages. These remote sensing derivate phenology dates were compared to National Agricultural Statistics Service (NASS) weekly crop progress reports to evaluate the capability of these datasets in temporal sensitive studies. The results illustrated that daily surface reflectance datasets were the most accurate source for time-sensitive studies. However, extra ancillary datum and appropriate denoising techniques should be applied to reconstruct the time series curve. Phenology matching process is a necessary step before detecting phenological information from remote sensing imagery for specific land cover type since same phenological stages of different crop types might have different counterparts on time series curve.

Introduction

Remote sensing based phenology detection techniques offer great opportunity for many climate change related researches, such as land surface vegetation dynamics [START_REF] Beck | Improved monitoring of vegetation dynamics at very high latitudes: a new method using MODIS NDVI[END_REF], agriculture management [START_REF] Garrigues | Modeling temporal changes in surface spatial heterogeneity over an agricultural site[END_REF], forestry monitoring [START_REF] Reference Ahl | Monitoring spring canopy phenology of a deciduous broadleaf forest using MODIS[END_REF] and global warming [START_REF] Parmesan | Influences of species, latitudes and methodologies on estimates of phenological response to global warming[END_REF]. Phenology is the study of the timing of cyclical events in nature such as springtime vegetation budburst or seasonal bird migrations [START_REF] Badeck | Responses of spring phenology to climate change[END_REF], while crop phenology refers to the developments, differentiation and initiation of organs of a crop (Hodges, 1991). Many researchers have studied phenological change using fieldcollected data, which is prohibitive for large scale monitoring due to the need for repeat data collection in order to monitor change from previous years (Thayn and Price, 2008). The situation has been totally changed after satellite sensor imagery, which record multi-scale (spatial or temporal) information of ground surface, was applied to monitor land vegetation dynamics.

A number of research works have been conducted to monitor "generic or mixed" vegetation phenological events [START_REF] Justice | Analysis of the phenology of global vegetation using meteorological satellite data[END_REF][START_REF] Reed | Measuring phenological variability from satellite imagery[END_REF][START_REF] Zhang | Monitoring vegetation phenology using MODIS[END_REF][START_REF] Fisher | Green leaf phenology at Landsat resolution: Scaling from the field to the satellite[END_REF]White and Nemani, 2006;[START_REF] Stöckli | Remote sensing data assimilation for a prognostic phenology mode[END_REF], and many efforts have been made to develop a general procedure for remote sensing based crop phenology detection. Sakamoto et al. [START_REF] Sakamoto | A crop phenology detection method using time-series MODIS data[END_REF] generated planted, heading and harvest dates of rice paddy in Japan using time series MODIS Enhanced Vegetation Index (EVI) datasets. Vina et al. [START_REF] Vina | Monitoring maize (Zea mays L.) phenology with remote sensing[END_REF] presented a visible atmospherically resistant index (VARI) and accumulated growing degree days (AGDD) to detect maize physiological transitions that undetectable from NDVI. Duchemin et al. [START_REF] Duchemin | Monitoring wheat phenology and irrigation in Central Morocco: On the use of relationships between evapotranspiration, crops coefficients, leaf area index and remotely-sensed vegetation indices[END_REF] monitored wheat phenology through the analysis of relationship among evapotranspiration, crop coefficients, leaf area index (LAI) and NDVI. Islam and Bala (Islam and Bala, 2008) made use of NDVI and LAI to evaluate potato phenological characteristics by modeling its essential phenological metrics. Each of these pioneering works has its special contribution to remotely-sensed crop phenology estimation research.

However, the majority of above-mentioned studies utilized multi-temporal NDVI/EVI/LAI values processed with the MVC technique [START_REF] Holben | Characteristics of maximum-value composite images from temporal AVHRR data[END_REF], which effectively removed the cloud-contaminated data points in original observations series. The maximum value of the composited period (usually range from 8 to 16 days) was assigned to this multiday period rather than the date when this value was captured (Thayn and Price, 2008). This processing step imbedded a temporal error of 8 (or 4) days in average and 15 (or 7) days maximum because a maximum value in 16 (or 8) days is selected and assigned to all of 16 (or 8) days. This could not be neglected in time-sensitive studies, especially for crop phenology detection. Therefore, it is essential to assess the effect of this composited process on crop phenology estimation accuracy and seek out an appropriate strategy for the remote sensing based crop phenology estimation.

The objectives of this paper are to: a) check the feasibility of employing daily surface reflectance datasets and its corresponding ancillary data to generate NDVI

Data Description

Study sites

In this study, the State of Iowa (range from 40°23′N to 43°30′N and 90°8 ′W to 96°38′W) in U.S. Corn Belt is selected for experiment because it is a predominant corn and soybean state and it is one of the most important corn and soybean production areas in United States. Iowa locates in the middle of American and its terrain is mainly dominated by plain. All these facts will help us to simplify the study and to keep focus on the method for phenological information extraction from the MODIS products. In this study, 50 sites (25 sites each for corn and soybean) within Iowa State, as shown in figure 1, were selected for method and data experiments.

Datasets description and preprocessing

Three MODIS products have been used in this experiment: Surface Reflectance Daily L2G Global 250m (MOD09GQ), Surface Reflectance 8-Day L3 Global 250m (MOD09Q1) and Vegetation Indices 16-Day L3 Global 250m (MOD13Q1). The MOD09GQ product provides an estimate of the surface spectral reflectance as it would be measured at ground level in the absence of atmospheric scattering or absorption. Low-level data can be corrected for atmospheric gases and aerosols, yielding a level-2 basis for several higher-order gridded level-2 (L2G) and level-3 products. Bands 1 (RED) and 2 (near infrared (NIR)) are provided at a 250-meter resolution in a daily gridded L2G product in the Sinusoidal projection. This product is meant to be used in conjunction with the Surface Reflectance Daily L2G Global 1km and 500m (MOD09GA) where important quality and viewing geometry information is stored. MOD09Q1 is a 8-day gridded level-3 product that stores an estimate of the surface spectral reflectance. Each MOD09Q1 pixel contains the best possible L2G observation during an 8-day period as selected on the basis of high observation coverage, low view angle, the absence of clouds or cloud shadow, and aerosol loading. MOD09GQ and MOD09Q1 are first calculated NDVI through the following formula: NDVI = (NIR-RED) / (NIR+RED). MOD13Q1, unlike MOD09GQ and MOD09Q1, is designed to provide consistent spatial and temporal comparisons of vegetation conditions. Blue, red, and near-infrared reflectance are used to determine the MODIS daily vegetation indices. The MODIS NDVI complements NOAA's Advanced Very High Resolution Radiometer (AVHRR) NDVI products and provides continuity for time series historical applications. Although MODIS also includes EVI that minimizes canopy background variations and maintains sensitivity over dense vegetation conditions, we only extract NDVI values for comparison purpose (USGS, 2009). All these datasets are captured between March 31 and December 11 from 2006 to 2008. This period includes 256 yeardays and totally covers the growing season for both corn and soybean in each year according to the statistical information published by United State Department of Agriculture (USDA). Many criterions, besides cloud status information recorded in quality flag, have been presented previously to remove noisy data points affected by thick clouds, gaseous and aerosol scattering as well as bidirectional reflectance distribution function (BRDF) distortion. For example, Sakamoto et al. [START_REF] Sakamoto | A crop phenology detection method using time-series MODIS data[END_REF] considered the pixels whose blue band reflectance was greater than 10% or whose sensor zenith angle was more than 32.25°should be removed as abnormal data. Points with a random NDVI increase greater than 0.4 during 20 days are also rejected as misleading data [START_REF] Chen | A simple method for reconstructing a high-quality NDVI time-series data set based on the Savitzky -Golay filter[END_REF]. We adopt these extra criterions to reserve the pixels in time-series of a test site for above-mentioned three datasets: a) the cloud and cloud shadow status should be "clear"; b) the blue band value should be lower than 0.1; c) the sensor zenith angle is smaller than 30°and d) Evaluation of temporal resolution effect in remote sensing based crop phenology detection studies 5

difference between survival data point and its Inverse Distance Weighted (IDW) estimation should be small than 0.1 [START_REF] Lu | Removal of Noise by Wavelet Method to Generate High Quality Temporal Data of Terrestrial MODIS Products[END_REF]. These standards could ensure that the survival sensor recorded reflectance values are not influenced by cloud or other noises. The National Agricultural Statistics Service (NASS), the USDA data collection arm,conducts hundreds of surveys every year and prepares weekly crop progress report containing crop planting and emergence percentage estimates. Crop progress data to be used for validation were acquired from the weekly 2006-2008 Crop Progress reports produced by the NASS for Iowa. (These reports are publicly available at http://www.nass.usda.gov/Statistics_by_State/Iowa/Publications/Crop_Progress_ &_Condition/index.asp.) According to these reports, there are nine districts in Iowa State: northwest (NW), north center (NC), northeast (NE), west center (WC), center (C), east center (EC), southwest (SW), south center (SC) and southeast (SE). Each of these districts has its own validation information (planted, emerged, silking, etc.) for different kinds of crops. Each stage was described as a period as its percentage covers from 0 to 100. Generally, the 50 percent corresponding Julian date of each stage could be selected as final validation information [START_REF] Wardlow | Using USDA crop progress data for the evaluation of greenup onset date calculated from MODIS 250-meter data[END_REF]. For corn, the report offered nine phenology stages: planted, emerged, tasseled, silking, milk, dough, dent, mature and harvest. For soybean, there are seven stages: planted, emerged, blooming, setting pods, leaf turning color, mature and harvest. All phenological stage definitions could be accessed at http://www.nass.usda.gov/Charts_and_Maps/Crop_Progress_&_Condition/Terms/ index.asp.

The Cropland Data Layer (CDL) program (http://www.nass.usda.gov/research/Cropland/SARS1a.htm), also offered by NASS, provides the information about what kind of crop has been planted on a specific field with coordinate information. USDA NASS online products are the Hu Zhao 1,2, Zhengwei Yang 3, Liping Di 2, Zhiyuan Pei1 main source of our experiment validation information.

Methodologies

Curve fitting strategy for daily data

After preprocessing steps, the time series are usually not equal-spaced or continuous since some abnormal pixels have been removed (as shown in Figure 2). Meanwhile, even those survival points are still contain some undefined noise (soil background, water content, moisture, saturation, pre-crop vegetation, etc.). Curve fitting is proved to be one of the most effective strategies to eliminate those noises [START_REF] Bradley | A curve fitting procedure to derive inter-annual phenologies from time series of noisy satellite NDVI data[END_REF].

Various techniques have been implemented to model the survival data points and to reconstruct the removed values, such as Best Index Slope Extraction (BISE) [START_REF] Viovy | The Best Index Slope Extraction (BISE): A method for reducing noise in NDVI time-series[END_REF], Fourier analysis [START_REF] Roerink | Reconstructing cloud free NDVI composites using Fourier analysis of time series[END_REF][START_REF] Geerken | Classifying rangeland vegetation type and coverage from NDVI time series using Fourier Filtered Cycle Similarity[END_REF], Savitzky-Golay filter and its improvement (Savitzky and Golay, 1964;[START_REF] Chen | A simple method for reconstructing a high-quality NDVI time-series data set based on the Savitzky -Golay filter[END_REF], Support Vector Machine [START_REF] Sun | Nonlinear Noise Reduction of Chaotic Time Series Based on Multi-dimensional Recurrent Least Squares Support Vector Machines[END_REF], asymmetric Gaussian function (Jönsson and Eklundh, 2002), mean-value iteration (MVI) (Ma and Veroustraete, 2006), wavelet transformation and its improvements [START_REF] Sakamoto | A crop phenology detection method using time-series MODIS data[END_REF][START_REF] Lu | Removal of Noise by Wavelet Method to Generate High Quality Temporal Data of Terrestrial MODIS Products[END_REF] and logistic function [START_REF] Zhang | Monitoring vegetation phenology using MODIS[END_REF]Jönsson and Eklundh, 2004;[START_REF] Beck | Improved monitoring of vegetation dynamics at very high latitudes: a new method using MODIS NDVI[END_REF][START_REF] Beck | A ground-validated NDVI dataset for monitoring vegetation dynamics and mapping phenology in Fennoscandia and the Kola peninsula[END_REF]. Hird and McDermid (Hird and McDermid, 2009) made a model-based empirical comparison of some frequently-used NDVI time series denoising techniques and concluded that double logistic function fitting strategy demonstrated the superiority over other methods. The formula could be expressed as equation ( 1), where wNDVI stands for minimum NDVI value (often arose in winter) of the final curve fitting time series. In this test, wNDVI represents bare soil or pre-crop vegetation NDVI value. This value is often represented as geometric mean of NDVI values in non-growing season (winter or bare soil background) [START_REF] Beck | Improved monitoring of vegetation dynamics at very high latitudes: a new method using MODIS NDVI[END_REF]. According to the NASS crop progress reports; there is no crop or plant available in the corn or soybean fields before April from 2006-2008. This means that we could set the non-growing season as the date before the 120th day of the year. mNDVI denotes the maximum NDVI and its difference with wNDVI represents the amplitude of NDVI time series; mS and S are variables relating to the inflection point in vegetative growing season, while mA and A are referred to senescence season.

Their relationship as shown in Figure 3 could be expressed [START_REF] Beck | A ground-validated NDVI dataset for monitoring vegetation dynamics and mapping phenology in Fennoscandia and the Kola peninsula[END_REF]. Least square method is commonly applied in function fitting to select the "best fit" parameters [START_REF] Steinier | Comments on Smoothing and differentiation of data by simplified least square procedure[END_REF][START_REF] Sun | Nonlinear Noise Reduction of Chaotic Time Series Based on Multi-dimensional Recurrent Least Squares Support Vector Machines[END_REF][START_REF] Markwardt | Non-linear Least Squares Fitting in IDL with MPFIT[END_REF] and often expressed as follows:

( ) ( ) 1 1 ( ) ( 1) 1 1 t mS t S mA t A NDVI wNDVI mNDVI wNDVI e e            (
  2 2 1 n i i ti i w NDVI NDVI       (2)
Where i w is the weight of ith NDVI value. When 2  reaches its minimum, the parameters are optimal. Beck et al. [START_REF] Beck | Improved monitoring of vegetation dynamics at very high latitudes: a new method using MODIS NDVI[END_REF] presented a two-step schema to determine weight, which was supposed to get the upper envelop of the survival data. First, all weights were set equally as 1 to obtain the reference curve line. In the second step, the weights were reassigned different values according to relative relationship between original values and the reference curve line: if the original value is above (below) the reference curve line, it gets a higher (lower) weight.

i NDVI is the ith originally calculated NDVI value, and ti NDVI is the corresponding estimated value. Since there are contaminated data in daily NDVI, we have made a little adjustment to this schema: set all contaminated points as zero and then applied Beck"s two-step weighting strategy. The weights of those contaminated data were set to zero throughout the whole weighting process. This ensures that those contaminated points will not affect the phenology estimation. To get more accurate phenology estimation, the composited products have been linearly interpolated to reconstruct contaminated daily data with ENVI Interactive Data Language (IDL) routine "interpol". However, all interpolated values are weighted as zero during the weighting process to keep the "original" survival Hu Zhao 1,2, Zhengwei Yang 3, Liping Di 2, Zhiyuan Pei1 points" affect on resultant NDVI curve.

Phenology stage estimation based on NDVI curve feature points

As shown in Figure 3, a NDVI curve could be characterized with seven curve feature points: G, S, P, M, Q, A, E. The inflection points S and A could be obtained from the maximum and minimum of the first derivate of estimated NDVI curve. The point M represents the date of the maximum NDVI value. If there exist more than one date, M will be the central point of these maximum-value dates. Points S, P, M, Q, and A divided the estimated NDVI curve into four parts. In each of these four parts, the maximum curvature could be obtained and marked as G, P, Q and E according to the following formula:

    3/ 2 2 1 t t NDVI NDVI      (3)
Julian date t is the independent variable. These seven temporal points will be used to match with actual phenological dates.

After calculating seven feature point dates for each test site, these dates should be matched with statistical stages to identify the relationship between curvederived dates and specific statistical phenological event dates (treated as actual phenology stages). Before matching, these remote sensing derived dates could be grouped as set G, set S, set P, set M and set Q. We also marked the statistical validation phenology dates (planted, emerged, silking, milk, mature and harvest, etc.) as phenology stage A, stage B, stage C, stage D, etc. By calculating rootmean-square error (RMSE) between each set and stage, an error matrix was produced through formula (4), in which εi,j was the RMSE between ith set , ik Set and jth phenology stage Stagej (includes m test sites, k was the serial number of test sites).

 

2 , 1 , m i k j k ij Set Stage m      (4)
All RMSE matrices for nine districts were first computed with test sites that falls into each district, and the state-level RMSE matrix was then obtained by weighting every district matrix. The weight of each district-level matrix is the proportion of district test site numbers to the total test sites. Figure 4 illustrates the whole process of final matching pair generation and error analysis.

Here are three basic principles to determine the intermediate match pairs for each dataset: 1) smaller RMSE, the minimum εi,j was selected to determine the first match pair, for example, Setp and Stageq; 2) time line consistency, the sets

Evaluation methods

We used two complementary approaches to evaluate the performance of each data set: phenology matching accuracy and phenology estimation error. First, we define three intermediate matching results and the final matching pair as: Imp01, Imp08, Imp16 and Fmp. They record the matching pairs as phenology "milk" to metric "M"", phenology "mature" to metric "A"", etc. We interpreted phenology matching accuracy (εi,j) as ratio between number of final matching pairs (Fmp) and the number of intersection of intermediate matching pairs (Imp01, Imp08 or Imp16) and the final matching pairs (Fmp). The matching accuracy reflects the degree of matching between intermediate matching results and the final matching pairs, εmi equal to 1 means the intermediate result is the same as the final result.

The bigger the εmi is, the worse they coincide with each other. The formula could be expressed as:

() () mi i N Fmp N Fmp Imp    (5)
where εmi stands for the phenology matching accuracy of datasets i, N is the function that extract the number of a set. Impi represents Imp01, Imp08 or Imp16 and symbol  means the intersection of the two datasets.

The phenology estimation error (εi,j) refers to the summation of final matching pair estimation errors, which is determined by the ratio of RMSE and the mean of the matching pairs" RMSE. For phenology estimation error εi,j, it happened between different datasets and could be listed as the ratio of final matching RMSE value RMSEij to the mean of its counterparts RMSEij (of different datasets). Formula 6 defines phenology estimation error.

m ij ei j ij i RMSE RMSE        (6)
The final error is the summation of matching error and estimation error. And it can be expressed as:

i mi ei ab         (7)
Where a and b are the weight for phenology matching error and phenology estimation error. Here we consider these two kinds of error have the same contribution to the final total error, which sets both a and b as 1.0.

Results and Discussion

According to the basic rules described in section 3.2, we got the final matching pairs for these three different temporal resolution datasets in 2006-2008 (figure 5 ). For corn, the final matching pairs in 2007 and 2008 are the same: field observation phenology stage "planted" corresponds to remote sensing derived first metric "G"; emerged to "S"; tasseled to "P"; milked to "M"; dough to "Q"; mature to "A" and harvest to "E". There is no metrics corresponded to phenology stage corn silking and dent in these two years. However, in 2006, the matching pairs seem slightly different: the stage emerged refers to metric "G" and there is no stage matching with metrics "S"; and for metric "A", the correspondent phenological stage is dent, not dough. Other matching pairs are the same with the following two years. According to the statistical observation information obtained Evaluation of temporal resolution effect in remote sensing based crop phenology detection studies 11 from NASS, the maize dough phenology stage is recorded 5 (South center of IOWA, 2007) to 15 (North center of IOWA, 2007) days before dent stage. When comparing to the RMSE error between remote sensing derived phenological dates and the validation information (in table 1 and2), the duration between dough and dent is somewhat negligible. This indicates that the matching results are not varying greatly across different years since phenology changes every year.

Comparing with vegetation types, comprehensive human being participation (seed selection, field management, irrigation status, etc.) is the most important factor of crop land cover development, which is sensitive to both natural and manmade microenvironment. For soybean, year 2006 and 2008 have the same matching pairs: planted to "G"; emerged to "S"; blooming to "P"; setting pods to "M"; leaves turning color to "Q" and harvest to "A". Remote sensing derived metric "A" refers to dropping leaf rather than harvest in 2007. The inter-annual matching error often happened at beginning or end of growing season. We believe that, besides the different cultivate conditions mentioned above, the limitation of NDVI at low vegetation density might be responsible for this experiment matching results. Another matter of concern was that, even using the same definition, stage "harvest" is matched with metric "E" for corn and "A" for soybean. This indicates that soybean was harvested when NDVI descended dramatically. This is consistent with the reality: soybean is harvestable after it dries completely (soybean has less than 15% moisture). At this stage, soybean plant turns completely brown. However, when corn reaches harvestable stage, its stem is usually still green. For corn and soybean, the 2007 state level RMSE results are listed in table 1 and 2 respectively, in which normal italic RMSE values (with or without underlines) are intermediate matching pairs for each datasets, while final matching results are highlighted by underlines. The RMSE range for corn is from 7.6 (between planted and "G", daily surface reflectance) to 25.6 (between emerged and "S", 8-day composited surface reflectance) days. The average difference between NASS observation and remote sensing results for these three datasets are 13.5, 14.4 and 13.9 days. For soybean, the RMSE ranges from 3.0 (between emerged and "S", 16 day composited vegetation product) to 23.9 (between turning color and "Q", 16 day composited vegetation product) days. And the average difference between statistical phenology stages and NDVI curve based estimates for these datasets are 14.4, Evaluation of temporal resolution effect in remote sensing based crop phenology detection studies 13

14.5 and 16.5 days. MODIS Daily surface reflectance improves the remote sensing estimation results in 1-2 days when compared with MODIS MVC products. The total errors for three datasets are displayed in figure 6. Generally the daily surface reflectance product performs better than other two datasets when being applied to crop phenology estimation studies in 2006-2008. The 8-day composited surface reflectance dataset works better than MODIS vegetation production (MOD13Q1) overall except a little irregularity for corn in 2007. Moreover, there is no obvious conclusion could be reached for crop-specific difference between corn and soybean.

Conclusion

Remote sensing based phenology detection is highly related to the temporal resolution of datasets. We have evaluated the effect of different datum which has different temporal resolution for two crop types: corn and soybean. USDA NASS crop progress reports have been used as ground truth data to assess the performance of different datum. Several conclusions could be drawn from above analysis: 1) Daily surface reflectance dataset could be employed to detect phenology information and had slight advantage to the composited data products, however, extra ancillary data and appropriate denoising and reconstruct techniques should be applied to eliminate the contaminated points. Double logistic function fitting method would be highly recommended; 2) The MVC process brought additional temporal error to final time series trend curve by removing contaminated points. This had a negative effect on time-sensitive remote sensing applications, for example, crop phenology estimation; 3) phenology matching is a necessary step when detecting specific stage for specific Hu Zhao 1,2, Zhengwei Yang 3, Liping Di 2, Zhiyuan Pei1 land cover type. Even some phenological stage have the same name, for example, emerged or harvest, they might corresponded to different segments of the final time series curve because of differences in cropping practice. Furthermore, phenology is a multi-factor dependent index to describe plant growth; NDVI has its own limitations to work as the only factor, other factors, like precipitation, temperature and soil moisture, et al. should be considered as part of the phenology model.

Fig. 1 .

 1 Fig. 1. Distribution of test sites in Iowa State, United States . Filled symbol pentagram (★) and triangle ( ▲ ) represented corn and soybean respectively. The background color gradient indicates the general condition of elevation.

Fig. 2 .

 2 Fig. 2. Time-series datasets of MOD09GQ (open circle 'o'), MOD09Q1 (dotted line) and MOD13Q1 (solid line) after preprocessing. Many MOD09GQ cloud-contaminated data points are removed and replaced with NAN so the time series looks like discrete.

1 )Fig. 3 .

 13 Fig. 3. Example of the double logistic function

Evaluation

  of temporal resolution effect in remote sensing based crop phenology detection studies 9 before (or after) Setp will only be matched with stages before (or after) Stageq, which is proposed to keep the consistency over time line. The final matching pair is generated by considering daily, 8-day and 16-day intermediate matching results, which should abide by the 3) majority rule: first we pick out the matching pairs in which these three datasets have three same intermediate results. Then we choose the matching pairs that two of these three datasets have same matching pairs and the other dataset has different result which will be re-decided same as the former two counterparts. If one stage refers to three different datasets, the result will be determined by the minimum RMSE value. The second rule, time line consistency, is simultaneously applied during this process.

Fig. 4 .

 4 Fig. 4. Flowchart of matching pairs determination and error analysis. It includes 5 steps, step 1: calculating district level RMSE matrix, step 2: integrating 9 district results into state level RMSE matrix; step 3: determining intermediate matching pairs for daily, 8-day and 16-day datasets; step 4: reaching a final matching pair for corn and soybean and step 5: conducting error analysis by using matching pair results in step 3 and 4.

Fig. 5 .

 5 Fig. 5. Final matching pairs for corn (left) and soybean (right) in 2006-2008.

Fig. 6 .

 6 Fig. 6. RMSE result for corn and soybean calculates form different datasets.

Table 1 .

 1 State level RMSE result for corn in 2007 (days)

		Planted Emerged tasseled Silking	milk	Dough	Dent Mature Harvest
	01G	7.6	11.3	64. 9	69.4	84.2	96.2	108. 6	127.4	171.0
	01S 26.4	16.6	39. 8	44.3	59.1	71.1	83. 5	102.4	144.4
	01P	53.1	43.3	15.3	19.3	33.3	45.1	57.3	76.1	118.0
	01M 78.0	68.2	15.4	12.3	11.4	21.0	32. 9	51.4	93.2
	01Q 99.3	89.4	34.5	30.1	16.3	9.9	15.5	31.4	72.4
	01A 133.3	123.4	68.0	63.6	48.9	37.1	25.3	11.5	39.3
	01E 127.9	125.7	94. 6	90.1	75.4	63.6	51.6	33. 7	22.8
	08G 13. 8	9. 6	55. 8	60.3	75.0	87.0	99.3	118.2	150.2
	08S 35.3	25.6	31.8	36.2	50.9	62.8	75.1	94.0	135.9
	08P	61.1	51.3	11.5	13.9	26.7	38.1	50.1	68.7	110. 5
	08M 85.2	75.4	21. 4	17.4	9.8	15.7	26.5	44.6	86.3
	08Q 103.7	93.7	38.7	34.3	20.2	10.6	10.6	26.5	67.7
	08A 135.6	125.7	70.2	65.7	51.0	39.1	27.0	10.5	36.2
	08E 140.8	134.5	95.1	90.6	75.9	64.0	51.9	33.7	19.7
	16G	8.2	14.1	68.7	73.2	88.1	100.1	112.5	131.4	152.0
	16S 22.0	12.8	44.4	48.9	63.7	75.8	88.1	107.0	149.0
	16P 46. 7	36.9	20.9	25.1	39.6	51.5	63.7	82.5	124.5
	16M 76.4	66.6	15.8	13.3	15.4	24.5	54.0	95. 6	35.8
	16Q 97.5	87.7	34.0	30.1	17.7	11.5	17.5	34.4	75.3
	16A 129.5	119.6	64. 2	59. 7	44.9	33.1	21.2	7.8	42.2
	16E 149.4	144.1	88. 6	84.1	69.3	57.3	45.2	27.0	21.2


Zhengguo Li
email: zhengguoli@caas.net.cn

Peng Yang
email: yangp@mail.caas.net.cn

Yanglin Wang
email: ylwang@urban.pku.edu.cn

Qinbo Zhou
email: zhouqb@mail.caas.net.cn

Huajun Tang
email: hjtang@mail.caas.net.cn

Hsiaofei Chang

Agricultural Landscape Dynamics and its Response in Seasonal Vegetation Activities in the Loess Plateau, Northern Shaanxi, China

Keywords: Agricultural Landscape, Landscape structure, Normalized Difference Vegetation Index (NDVI), The Loess Plateau of China

The ecological and environmental conditions in semiarid areas are closely linked to landscape dynamics. This study examined the seasonal vegetation activities of landscape classes and dynamics in the Loess Plateau in North Shaanxi Province, China. Landscape dynamics were studied by comparing Landsat images in 1987 and 2007, and classified as landscape transitions, landscape variations and unchanged landscape. The characteristics of seasonal vegetation activities, such as dates for onset-of-greenness, peaks, and onset-of-senescence for landscapes were determined using Moderate Resolution Imaging Spectroradiometer (MODIS) 10-day vegetation index composite products (2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009). The landscape dynamics showed that some sloping croplands were converted to woodlands and grasslands. Vegetation activities responded diversely and seasonally to landscape classes and dynamics. Temporal analysis of the Normalized Difference Vegetation Index (NDVI) showed that the time of onset-of-greenness and onset-of-senescence were similar for all major landscape classes, while the time of peak vegetation activity was longer for tree than grass species. The analysis of seasonal vegetation activities and landscape dynamics indicated that seasonal vegetation activities were not only relative to landscape classes, but also affected by landscape dynamics.

Introduction

Serious soil erosion in Loess Plateau area of China and the increase in human disturbances in ecosystem have resulted in a range of ecological and land use problems, such as the depletion of land resources and the degradation of ecoenvironments. An urgent need for agricultural landscape and ecological planning appears to ensure the sustainable development of the Loess Plateau area (Fu et al., plant forms [START_REF] Fu | Agricultural landscape spatial pattern analysis in the semiarid area of the Loess Plateau, China[END_REF][START_REF] Fu | Temporal change in land use and its relationship to slope degree and soil type in a small catchment on the Loess Plateau of China[END_REF]. Within this system, the time series of NDVI can characterize the temporal extent of the growing season and productive potential of agricultural landscape at the regional scale.

Capability to assess the relationship between landscape dynamics and vegetation seasonal changes is required for recognizing the critical factors and modelling the ecological condition changes in the Loess Plateau area [START_REF] Xiao | Ecological Principles for landscape classification and evaluation[END_REF][START_REF] Fu | Agricultural landscape spatial pattern analysis in the semiarid area of the Loess Plateau, China[END_REF]. The vegetation response to landscape dynamics can reveal any benefits or damage to eco-functions of the regional ecosystem [START_REF] Wang | The spatial structure of landscape eco-systems: concept, indices and case studies[END_REF].

The landscape dynamics and its response in seasonal vegetation activities in the Loess Plateau area, China were analysed for a period of 20 years. The specific objectives were: (1) to investigate the dynamic changes in the agricultural landscape and identify the main driving factors; (2) to determine the periods for onset-ofgreenness, peaks and onset-of-senescence in each agricultural landscape, [START_REF]Institute of Natural Resources and Regional Planning[END_REF] to quantitatively and qualitatively assess the sensitivity of NDVI in discriminating landscape dynamics at the regional scale.

Study area and materials

Study area

The study area was the northern part of Shaanxi Province, China between lat 35°2 1′N-39°34′N and long 107°28′E-111°15′E, in the middle part of the Loess Plateau. The area is 80,606 km 2 and the main rivers include the Kuye, Wuding, Qingqian, Yan, Fenchuan Rivers and part of the Luo River (Fig. 1). The terrain has significant topographic variability with typical Loess hills and gully slopes with elevations of 400-1900 m. The area has a semi-arid continental monsoon climate, with cold winters and warm summers. During the past decades, agriculture has been the primary economic activity and income source in the area, and the most important land use change was the clearing of forestland and grassland for agriculture [START_REF] Chen | Land use change in a small catchment of northern Loess Plateau, China[END_REF]. The hilly topography, intensive precipitation and especially long-term extensive human activity (i.e. removal of the natural vegetation and farmingaccelerated deterioration) have caused serious soil erosion [START_REF] Fu | Agricultural landscape spatial pattern analysis in the semiarid area of the Loess Plateau, China[END_REF][START_REF] Fu | Temporal change in land use and its relationship to slope degree and soil type in a small catchment on the Loess Plateau of China[END_REF]. The present land cover and vegetation types include cropland, shrub, woodland (mainly in southern mountainous areas) and grassland. The main crops are potatoes (Solanum tuberosum), beans (Phaseolus vulgaris), maize (Zea mays) and millet (Panicum miliaceum). The woodland is dominated by introduced vegetation, mainly locust trees (Robinia pseudoacacia). In this area, little leaf peashrub (Caragana microphylla) is the most important shrub species. The grassland is mainly covered by annuals such as wheatgrass (Agropyron cristatum Gaertn.), sweet wormwood (Artemisia annua), annual fleabane (Erigeron annuus) and sandy needlegrass (Stipa glareosa).

Data used

The landscape class data were derived from a series of Landsat-5 TM images acquired in 1987 and 2007. In total there were 12 images from each of the two time periods with the orbit sequence numbers 126-32, 34, 35 and 36; 127-33, 34, 35 and 36; and 128-34 and 35. MODIS data were also used. MODIS 10-day vegetation index (VI) composite products in the period 2002-2009 were applied in the analyses. Due to limited availability of validated Landsat data coinciding with MODIS collects, the analysis was limited to 2007-2009, and to the major landscape classes( Table 1).

Methods

Data pre-processing

Data pre-processing used remotely-sensed image processing software packages-ENVI 4.0, developed by RSI Company, and the geographical information system software ArcGIS 9.0, by ESRI Company. Prior to analysis, following atmospheric radiance calibration and initial geometric rectification of all the TM images using ENVI, a number of reference points were selected from a 1:50,000 scale topographic map. The TM images were rectified to a Gauss Kruger projection with a pixel resolution of 30 m × 30 m, by using nearest neighbour rules provided by ENVI, to ensure that the error was controlled in less than one-half pixel root mean square error. The MODIS products were also rectified to the Gauss Kruger projection by the above method, and resampled to 500 m × 500 m pixel resolution. Finally, the data was masked so that the analysis was limited to the administrative area of the northern Shaanxi Province.

The MODIS VI products included maximum value composites of the NDVI and enhanced VI, as well as the corresponding red and near-infrared reflectance (band 1 and band 2) and quality assessment (QA) flags [START_REF] Huete | Overview of the radiometric and biophysical performance of the MODIS vegetation indices[END_REF]. Only NDVI values from the MODIS composite products were used, along with the VI Usefulness Index in the QA data set. Reflectance values of NDVI for dates with a Usefulness Index value lower than -good quality‖ were replaced by linearly interpolated values from the two closest dates with good, high or perfect quality. Finally, an iterative median smoother was used to retain peaks and eliminate valleys of the temporal NDVI curve.

Classification of the 1987 and 2007 Landsat TM images

Landscape classification and mapping with remotely sensed data is a common technique, mainly based on the ability to reflect land surface and human activity information [START_REF] Ichoku | A review of mixture modeling techniques for sub pixel land cover estimation[END_REF][START_REF] Metternicht | Estimating erosion surface features by linear mixture modeling[END_REF]. We used supervised maximum likelihood classification methods. The accuracy of the 1987 and 2007 images was assessed using ground data points not used in the classification process. The equal control point methods in ENVI were used with at least 30 points for each class. The overall accuracy assessment was checked separately for each image and accepted if the accuracy was > 80% (determined with 580 random checkpoints). After accuracy assessment, all images were clumped and vectorized using ENVI. These coverages were pre-processed to eliminate areas < 0.27 ha (corresponding to 3 pixel × 1 pixel) using ArcGIS for faster spatial analysis.

Landsat TM images (1987) were successfully classified for 11 landscape classes; settlements and sparse woodland were classified with a lower accuracy than other classes. However, this was acceptable as the overall classification accuracy is much higher (87%) with the Kappa statistics value of 0.92. The classification of images (2007) was also acceptable due to a higher overall classification accuracy (89%) and Kappa values of 0.93. Three main agricultural landscape classes were identified: cropland, woodland and grassland. A further subclassification of agricultural landscape was defined using the percentage of vegetation canopy and terrain characteristics (Table 2).

Classifications of landscape dynamics

According to the intensity or range of landscape changes, landscape dynamics were simplified to three typical landscape changes. First, the original landscape was replaced by a different landscape, causing transition of landscape classes. Second, the landscape class was unchanged while landscape subclass changed with varying landscape characteristics. Third, the landscape was unchanged because neither landscape class nor subclass changed.

The complexity of landscape dynamics was based on the landscape classification results in 1987 and 2007. The landscape classification scheme was defined as unchanged landscape, landscape variation and landscape transition. In detail; unchanged landscape meant the class and subclass of landscape remained as original; transition meant a shift of landscape class; and variation meant landscape class changed but subclass remained unchanged(Table 3).

Table 3 Classifications of landscape dynamics

Classes

Cropland Landscape Woodland Landscape Grassland Landscape

Subclasses Cropland-I Cropland-T Cropland-S Woodland-D Woodland-M Woodland-S Grassland-D Grassland-M Grassland-S Cropland Landscape Cropland-I ☆ Cropland-T ◇ ☆ Cropland-S ◇ ◇ ☆ Woodland Landscape Woodland-D △ △ △ ☆ Woodland-M △ △ △ ◇ ☆ Woodland-S △ △ △ ◇ ◇ ☆ Grassland Landscape Grassland-D △ △ △ △ △ △ ☆ Grassland-M △ △ △ △ △ △ ◇ ☆ Grassland-S △ △ △ △ △ △ ◇ ◇ ☆
Notes: unchanged landscape, ◇ landscape variation, △landscape transition.

Characteristics of seasonal vegetation activities

A preliminary study on individual years indicated that although the relative magnitudes of the NDVI differed from year to year, partly due to different view angles and atmospheric conditions, the relative time periods for major vegetation activities could still be inferred [START_REF] Senay | Combining AVHRR-NDVI and landuse data to describe temporal and spatial dynamics of vegetation[END_REF]. It was anticipated that averaging of many years of data could eliminate sudden and unexplainable irregularities in the NDVI data [START_REF] Michael | Estimating spatio-temporal patterns of agricultural productivity in fragmented landscapes using AVHRR NDVI time series[END_REF]. Thus, the source data were from the 3-year average 10-day maximum value MODIS-NDVI composites (2007-2009).

The key to defining seasonal characteristics is to identify the onset of the growing season. The method used was visual inspection of the seasonal NDVI curves accompanied by an evaluation of the difference between adjacent NDVI values [START_REF] Senay | Combining AVHRR-NDVI and landuse data to describe temporal and spatial dynamics of vegetation[END_REF]. Once the onset of the growing season was identified, additional seasonal characteristics were derived (Table 2). There were more than one local maximum NDVI values; two local maxima occurred, in the beginning and later part of the growing season. There were two families of seasonal characteristics: temporal and NDVI value. Temporal characteristics are time of onset of growing season (OnT), time of end of growing season (EndT), time of the first maximum NDVI (MaxT-1), and time of the second maximum NDVI (MaxT-2). NDVI characteristics are NDVI value at onset of growing season (OnV), NDVI value at end of growing season (EndV), the first maximum NDVI value (MaxV-1), and the second maximum NDVI value (MaxV-2).

Results

Changes of the landscape pattern between 1987 and 2004

According to the classified Landsat images between 1987 and 2007, landscape changes have occurred in composition and structure of the Loess Plateau area (Table 5). Cropland was the dominant landscape class in this area. The proportion of cropland landscape decreased from 38 to 36.6%; especially, there was a net decline of 1128.5 km 2 in terrace cropland during 1987-2007. Grassland was also a significant landscape class, its area increased from 46.1 to 47.3%, with a net increase of 806.1 km 2 in sparse grassland. Overall, a large amount of cropland was transformed into woodlands and grasslands. Hence, over 1987-2004, there was a decrease in cropland and an increase in grassland and woodland.

Landscape dynamics and landscape classes

Landscape dynamics among the major landscape classes, included landscape transition, variation and unchanged, were determined based on the classification of Landsat images. Certain landscape transitions or variations were observed and evaluated for each class (Table 6). For example, nearly 80.9 km 2 (25.1%) of irrigated cropland was converted into other landscape classes; however, there were no landscape variations. A large amount of terrace cropland (882.2 km 2 ) was transformed into grassland or woodland landscape, and 318.6 km 2 of terrace cropland developed into irrigated cropland or degraded to sloping cropland. Of sloping cropland, a net area of 510.7 km 2 shifted into grassland or woodland, and 116.1 km 2 upgraded into irrigated or terrace cropland.

Of woodland landscape, the unchanged dense woodland was 1954.7 km 2 (97.0%), only 30.2 km 2 was converted into other landscape classes, and a similar area of dense woodland degraded into medium or low coverage woodland. Most moderate dense woodland remained unchanged (95.6%) and 1.0% changed into other woodlands, while 243.9 km 2 (3.4%) shifted into other landscape classes. Of sparse woodland, there was an area of 123.8 km 2 in landscape transitions and 15.5 km 2 in variations.

Of grassland landscape, 2011 km 2 (89.1%) of dense grassland remained unchanged, 243.8 km 2 (10.8%) was converted to cropland and 2.3 km 2 (0.1%) was degraded. A large amount of moderate dense grassland 3977.9 km2 (17.5%) was heavily degraded into low coverage, and 636.5 km 2 (2.8%) changed into other landscapes. As for sparse grassland, an area of 1411.9 km 2 was in landscape transitions and 3700.1 km 2 in variations.

Landscape dynamics and spatial distribution

The analysis of landscape composition in each watershed showed that the proportion of cropland landscape decreased gradually from north to south in the Loess Plateau area, while woodland and grassland obviously increased. Further landscape dynamics in watersheds were subsequently evaluated.

The landscape classes in the north, where the classes were cropland and grassland, such as in the watersheds of Gushan Chuan River, the Kuye River and the middle reaches of the Wuding River, generally did not change (Table 7). In the middle area, the variations within landscape classes were significant; in Toudao Chuan River (36.9%), the upper (18.2%) and middle (19.5%) reaches of the Wuding River, the upper (17.8%) reaches of the Yan River and the upper (17.1%) reaches of the Luo River, where the dominant landscape classes were dense and moderate dense grassland. Landscape transitions occurred mainly in the northwest, such as the continental watersheds (19.2%), Tuwei (15.3%) and Yuxi Rivers (13.3%), where the landscape classes were mainly sparse grassland and sloping cropland. In the south, where woodland and grassland were the major landscape classes, almost no landscape transition took place; such as the watersheds of the upper reaches of the Luo (0.8%), Ju (1.1%), Yunyan (1.7%), Shiwang (2.8%) and Hulu (3.2%) Rivers.

Seasonal vegetation activities of landscape subclasses

Based on temporal plots of the 3-year average NDVI, the seasonal vegetation activities (onset, peak and senescence) of the various landscape subclasses were determined. For each landscape subclass, onset-of-greenness usually occurred about 15 April (Table 8). The peak-of-greenness occurred at different periods for different subclasses. The irrigated cropland was the first cropland to peak at the end of May; this was expected since this subclass was dominated by winter wheat. MaxT-1 occurred on approximately 15 July for terrace and sloping cropland, sparse woodland, moderate dense and sparse grassland; for dense woodland and grassland, moderately dense woodland, it occurred around 15 May. MaxT-2 was about 5 September for all cropland, sparse woodland, moderate dense and sparse grassland; and on about 25 September for dense and moderate dense woodland and dense grassland. Onset-ofsenescence was about 25 October for all the subclasses.

Landscape subclasses are compared in Table 8. All landscape subclasses in the Loess Plateau area had a high degree of variability in NDVI. Average NDVI values ranged from 6.1 for terrace cropland to 18.3 for dense woodland at the OnT. Furthermore, NDVI ranged from 33.7 and 38.3 for sparse grassland to 62.2 and 61.2 for dense woodland at the MaxT-1 and MaxT-2, respectively. Similarly, the NDVI values varied from 30.0 to 57.3 at the EndT.

Seasonal vegetation activity of landscape dynamics

The temporal plots of NDVI showed that vegetation activities responded diversely and seasonally to landscape classes and dynamics. As for unchanged landscape, landscape subclasses in the Loess Plateau area could be categorized into two general groups of seasonal vegetation activities. The first group was irrigated cropland (Fig. 2a), dense and moderate dense woodland (Fig. 2b) and moderate dense grassland (Fig. 2c), which had an extended period between MaxT-1 and MaxT-2. The group could be influenced by the physiology of trees. The dense grassland and irrigated cropland subclasses also contained significant tree species inclusions (i.e. the MODIS pixels in a predominantly pasture or crop area might include signals from nearby forest). The second group was generally dominated by grass species. Sloping and terrace cropland (Fig. 2a), moderate dense and sparse grassland (Fig. 2c) had a relatively short period between MaxT-1 and MaxT-2. This could be due to these subclasses being mainly distributed in the northwest of the Loess Plateau area where rainfall is relatively low with a semi-arid climate. and moderate dense woodland landscape increased. The quantitative evidence showed that some sloping croplands were converted to woodlands and grasslands in the past decades. Presumably, this is due to structural adjustment of the agricultural sector, the initiation of Chinese Family Planning Policy, the increase in regional income in recent years and implementation of large-scale eco-environment conservation programmes (e.g. the Grain-for-Green Programme in 1999). From an ecological viewpoint, the present landscape structure and pattern is more sustainable than the former.

In semiarid and arid areas, vegetation has significant ecological functions, such as soil and water conservation and providing habitats for all life forms. Vegetation activities (onset, peaks and senescence) could directly reflect soil resistance to erosion and the temporal pattern of water use. Most of all, the identification of two growing season peaks could be important to identify water use patterns of different agricultural landscape types. The period between the second peak and senescence could signify a gradual decrease in water use, while the period after senescence a rapid decline in evapotranspiration [START_REF] Senay | Combining AVHRR-NDVI and landuse data to describe temporal and spatial dynamics of vegetation[END_REF].

The findings suggest that the spatial correspondence of MODIS and Landsat data sets was robust, and the combination of the two data sets to characterize spatial heterogeneity and temporal dynamics of vegetation was justified. This research has shown that MODIS NDVI data has the potential to determine periods of important vegetative activities and relative ‗greenness' magnitudes that can influence major energy and mass flux processes at regional scales.

Conclusion

On the regional scale, landscape compositions in the northern watersheds were mainly grassland and mixed cropland-grassland landscape, while in the southern watershed, landscape compositions mainly included mixed cropland-woodland and mixed woodland-grassland. Correspondingly, the vegetation cover ranged from low to high, which reflected the spatial pattern of vegetation cover resulting from natural forces and human activities on the macro regional scale.

Using a 3-year average of MODIS NDVI data, the temporal analysis of NDVI values showed that the onset-of-greenness (middle April) and onset-of-senescence (late October) were at similar times for all major landscape classes. The main difference among the classes was the NDVI magnitude in the period between onsetof-greenness and onset-of-senescence. Generally, the length of time for the peak vegetation activity (i.e. MaxT-1 to MaxT-2) was longer for tree than grass species.

The ranking of landscape subclasses in seasonal vegetation activities indicated that moderate dense grassland had landscape subclasses in two general groups. It is important to note that the moderate dense grassland was by far the dominant landscape in almost all watersheds. Landscape subclasses with seasonal average NDVI values less than moderate dense grassland included: all croplands, sparse woodland and sparse grassland; landscape subclasses with NDVI values larger than included dense and moderately dense woodland and dense grassland.

Taking landscape dynamics into account, the peak vegetation activity in unchanged landscapes was from high to low in the following sequence: dense and moderate dense woodland > dense grassland > sparse woodland and irrigated cropland > moderate dense grassland > terrace cropland > sloping cropland and sparse grassland. For landscape variations, the sequence was: dense and moderate dense woodland and irrigated cropland > dense grassland, sparse woodland and terrace cropland > moderate dense grassland > sparse grassland and sloping cropland. For landscape transitions, the sequence was: irrigated cropland > dense grassland > dense woodland, terrace cropland and moderate dense grassland > sloping cropland > moderate dense, sparse woodland and sparse grassland.

By analysing the relationships among the seasonal vegetation activities, landscape types and dynamics, it was evident that seasonal vegetation activities were not only relative to landscape classes, but also affected by landscape dynamics (e.g. landscape transitions and variations). Therefore, more attention should be paid to landscape spatial structure and temporal dynamics for the ecological recovery and construction in the Loess Plateau of North Shaanxi.

Fig. 1 .

 1 Fig.1. Geographical location of study area and delineation of watersheds Fig.1. Geographical location of study area and delineation of watersheds

Table 1 .

 1 Data source

	Data source		Date/resolution
		126-33/34/35/36	1987.07.29, 2007.06.02, with a resolution of 30m
	TM images	127-33/34/35	1987.06.02, 2007.06.25, with a resolution of 30m
		128-34/35	1987.05.24, 2007.08.03, with a resolution of 30m
	MODIS images	VI products (MOD13)	10-Day composite from 2007 to 2009, with a resolution of 0.5km
	Topographic map of Shaanxi province (1:50000)	Published by National Geomatics Center of China
			Made in 2000 by Institute of Geographical Sciences and Natural
	Land use map of Shaanxi province
			Resources Research, Chinese Academy of Sciences

Table 2

 2 Descriptions of landscape classes

	Landscape
		Landscape subclasses	Description
	classes
		irrigated cropland (Cropland-I)	cropland less than 15° (AD15)
			slope cropland more than 15° (AD15)
	Cropland landscape	terrace cropland (Cropland-T)
			and less than 25°(AD25)
		sloping cropland (Cropland-S)	slope cropland more than 25° (AD25)
		dense grassland (Grassland-D)	75-100% grass cover
	Grassland landscape	moderate dense grassland (Grassland-M)	50-75% grass cover
		sparse grassland (Grassland-S)	< 50% grass cover
		dense woodland (Woodland-D)	70-100% canopy cover
	Woodland landscape	moderate dense woodland (Woodland-M)	50-75% canopy cover
		sparse woodland (Woodland-S)	< 50% canopy cover
	Residential landscape	settlements	residential areas
	Water landscape	wetland	natural and artificial lake, river

Table 4

 4 Description of NDVI metrics calculated for the period 2007-2009

	Abbreviation Definition	Metric
		The period after which there was an increase of NDVI	Starting date of NDVI high period
	OnT
		values from a seasonally low value	(onset-of-greenness)
	OnV	Value of NDVI at forwards intersection	NDVI at start of high period
	MaxT-1	Time of the first maximum raw corrected NDVI	Date of the first maximum NDVI
	MaxV-1	The first Maximum value of NDVI	The first maximum of NDVI
	MaxT-2	Time of the second maximum raw corrected NDVI	Date of the second maximum NDVI
	MaxV-2	The second Maximum value of NDVI	The second maximum of NDVI
	EndT	the period after which a downward trend in NDVI began End date of NDVI high period
			NDVI at end of high period
	EndV	Value of NDVI at backwards intersection
			(onset-of-senescence)

Table 5

 5 Landscape structure compositions in study area

			Area (km 2 )			Ratio (%)
	Landscape type
		1987	2004	Difference	1987	2004	Difference
	Cropland-I	322.4	403.0	80.6	0.4	0.5	0.1
	Cropland-T	24504.2	23375.7	-1128.5	30.4	29.0	-1.4
	Cropland-S	5803.6	5723.0	-80.6	7.2	7.1	-0.2
	Woodland-D	2015.2	2095.8	80.6	2.5	2.6	0.1
	Woodland-M	7173.9	7335.1	161.2	8.9	9.1	0.2
	Woodland-S	2579.4	2498.8	-80.6	3.2	3.1	-0.1
	Grassland-D	2257.0	2095.8	-161.2	2.8	2.6	-0.2
	Grassland-M	22730.9	23053.3	322.4	28.2	28.6	0.4
	Grassland-S	12171.5	12977.6	806.1	15.1	16.1	1.0
	settlements	322.4	806.1	483.6	0.4	1.0	0.6
	wetland	725.5	241.8	-483.6	0.9	0.3	-0.6

Table 6

 6 Landscape dynamics of landscape types

			Area (km 2 )			Ratio (%)
	Landscape type	Unchanged	landscape	Landscape	Unchanged	landscape	Landscape
		landscape	variation	transition	landscape	variation	transition
	Cropland-I	241.5	0.0	80.9	74.9	0.0	25.1
	Cropland-T	23303.5	318.6	882.2	95.1	1.3	3.6
	Cropland-S	5176.8	116.1	510.7	89.2	2.0	8.8
	Woodland-D	1954.7	30.2	30.2	97.0	1.5	1.5
	Woodland-M	6858.2	71.7	243.9	95.6	1.0	3.4
	Woodland-S	2440.1	15.5	123.8	94.6	0.6	4.8
	Grassland-D	2011.0	2.3	243.8	89.1	0.1	10.8
	Grassland-M	18116.5	3977.9	636.5	79.7	17.5	2.8
	Grassland-S	7059.5	3700.1	1411.9	58.0	30.4	11.6

Table 7

 7 Composition of landscape dynamics in watersheds

			Ratio (%)
	Watersheds	Unchanged	landscape	Landscape
		landscape	variation	transition
	Gushan Chuan	98.8	0.2	1.0
	Kuye River	86.9	7.4	5.7
	Tuwei River	78.5	6.2	15.3
	Yuxi River	82.1	4.6	13.3
	Wuding River (upper)	72.4	18.2	9.4
	Wuding River (middle )	79.4	19.4	1.2
	Wuding River (lower)	94.3	1.5	4.2
	Jialu RIver	83.9	13.1	3.0
	Dali River	83.2	12.6	4.2
	Toudao Chuan	60.6	36.9	2.5
	Qingqian River	84.4	11.6	4.0
	Yan River (upper)	81.1	17.8	1.1
	Yan River (lower)	86.8	9.3	3.9
	Luo River (upper)	82.1	17.1	0.8
	Luo River (middle)	89.1	4.1	6.8
	Yunyan River	98.2	0.1	1.7
	Hulu River	92.4	4.4	3.2
	Shiwang River	96.2	1.0	2.8
	Ju River	97.7	1.2	1.1
	Continental watersheds	76.7	4.1	19.2

Table 8

 8 Summary of vegetation activity based on MODIS data(2007)(2008)(2009)

	Landscape class	OnT	OnV	MaxT-1	MaxV-1 MaxT-2	MaxV-2	EndT	EndV
	Cropland-I	4/15	13.5	5/25	41.2	9/5	55.3	10/25	48.3
	Cropland-T	4/15	6.1	7/15	34.5	9/5	40.3	10/25	30.9
	Cropland-S	4/15	8.6	7/15	39.6	9/5	43.8	10/25	33.8
	Woodland-D	4/15	18.3	5/15	62.2	9/25	61.2	10/25	57.3
	Woodland-M	4/15	14.2	5/15	47.5	9/25	53.1	10/25	49.7
	Woodland-S	4/15	11.9	7/15	50.6	9/5	50.5	10/25	44.2
	Grassland-D	4/15	16.3	5/15	58.6	9/25	59.8	10/25	56.6
	Grassland-M	4/15	8.6	7/15	39.9	9/5	43.5	10/25	36.0
	Grassland-S	4/15	7.5	7/15	33.7	9/5	38.3	10/25	30.0

Acknowledgements

The initial and ongoing research support was financially supported by National Basic Program of China (973 Program, 2010CB951502), by the National Natural Science Foundation of China (No. 40930101, No. 41001381 and 41001246), and by the Ministry of Finance of China through Non-profit National Research Institute (IARRP-2011-015). All persons and institutes who kindly provide their data available for this analysis are greatly appreciated.

Similar rules could be found in landscape variations. Differing to other croplands, the seasonal vegetation activities in irrigated cropland were characterized by higher maximum NDVI values at MaxT-1 and MaxT-2, shorter renewal period between OnT and MaxT-1, much more extended growing period between MaxT-1 and MaxT-2 (Fig. 3a). Similarly, the dense and moderate dense woodlands (Fig. 3b) and dense grasslands (Fig. 3c) had almost the same characteristics as each other of vegetation seasonal activities, while the sparse woodlands (Fig. 3b), moderate dense and sparse grasslands (Fig. 3c) had the opposite features with lower maximum NDVI values at MaxT-1 and MaxT-2, much longer renewal period between OnT and MaxT-1, shorter growing period between MaxT-1 and MaxT-2 Compared with irrigated cropland in landscape variations, there was a similar temporal pattern of seasonal vegetation activities in landscape transitions (Fig. 3d). However, the terrace and sloping croplands had lower maximum values of NDVI, longer renewal periods and shorter growing periods (Fig. 3d). Especially, MaxT-1 for these two subclasses was two months later than for irrigated cropland. At the same time, all woodlands converted from other landscape classes had much lower maximum NDVI values and shorter growing seasons than other woodlands (Fig. 3e). All dense grasslands in various landscape dynamics had high maximum values of NDVI and extended growing seasons (Fig. 3f), while other grasslands in landscape transitions and variations had the low maximum and short growing season

Discussion

The dynamics of landscape classes were identified by repeated satellite images coupled with GIS analyses of 1987 and 2007. The proportion of terrace and sloping cropland decreased. The proportion of moderate dense and sparse grassland, the dense


Zhihua Diao
email: diaozhua@163.com

Anping Zheng

Yuanyuan Wu
email: yuanyuan539@yahoo.com.cn

Shape Feature Extraction of Wheat Leaf Disease Based on Invariant Moment Theory

Keywords: feature extraction, shape feature, invariant moment, wheat disease

ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

With the development of technology and the improvement of computer image processing technology, image processing and pattern recognition technology can be used as an effective way for the diagnosis and identification of wheat diseases. Feature extraction is a key research topic in image processing and pattern recognition, and broadly refers to a transformation which transforms the samples in high dimensional space to low-dimensional space through the map or transformation to achieve the purpose of dimensionality reduction [START_REF] Wang | A Survey of Feature Selection[END_REF]. In wheat leaf disease aspects, shape feature is an important feature to describe the image. And as a local feature, it describes the geometric properties of the local area. As the target shape of images obtained from different vision may be involved in a large difference, translation, scale and rotation transformation invariance problems need to be solved for accurate shape matching.

In recent years, many researchers have carried out the image shape feature extraction research and achieved certain results in agricultural field. Chen Jiajuan et al [START_REF] Chen | Automatic Measurement of Danger Degree of Cotton Insect Pests Using Computer Vision[END_REF] extracted some shape features such as the length of the chain, chain code area and shape factor, tortuosity, irregularity and so on based on the holes and incomplete edges in cotton leaves to determine the extent of the victims of cotton pests. Xie Guojun et al [START_REF] Xie | Method for fruit shape feature acquisition based on multidirectional vision[END_REF] proposed a fruit shape feature extraction algorithm based on multidimensional vision system for fruit intelligent selection and grading. Cai Jianrong et al [START_REF] Cai | Shape feature extraction of on-tree citrus based on genetic algorithms[END_REF] obtained the best-fitting ellipse parameters to describe orange shape using genetic algorithm. Wang Xiaofeng et al [START_REF] Wang | Feature Extraction and Recognition for Leaf Images[END_REF] introduced a method to identify the leaf based on shape features of leaf images. Zhao Jinhui et al [START_REF] Zhao | Image Segmentation Method for Sugarcane Diseases Based on Color and Shape Features[END_REF] extracted the lesion area width, squareness and roundness for the separation of disease lesion according to sugarcane seedling image characteristics of red rot and ring spot diseases. Yang Zengfu et al [START_REF] Yang | Recognition of the Tea Sprout Based on Color and Shape Features[END_REF] detected the edge of tea buds according to the shape characteristics of tea buds, and which was used for tea buds identification. However, these studies did not consider the basic translation, scale and rotation transformation invariance problems, and had some flaws.

Invariant moment theory was introduced for shape feature extraction of wheat leaf disease image in this paper, and seven Hu Invariant moments were selected as shape parameters. Then shape features of wheat leaf disease image were extracted. A feature subset including twelve shape features was received, and was applied to wheat leaf disease identification system which acquires relatively high recognition efficiency.

Invariant Moment

Invariant moment is a statistical image feature which meets the translation, rotation and scale invariance, and has been widely used in image recognition. There are many invariant moments Used to describe shape which are arithmetic invariant moment, Legendre moment and Zerrnike moments, etc. [START_REF] Ding | Research on moment invariants algorithm[END_REF].

For the integral area of S , if we have a given two-dimensional continuous function ( , ) f x y , the () pq  order moment is defined as follows, ( , )

, , 0,1, 2, pq pq m

x y f x y dxdy p q

  . ( 1
)
If the digital image ( , ) f x y is a piecewise continuous function, and has non-zero value in a limited part of the xy plane, then we can prove the existence of various order moments. And the moment sequence pq m could be uniquely determined by the image ( , ) f x y . At the same time the moment sequence pq m also uniquely determine the image ( , ) f x y . The central moments can be expressed as follows,

( ) ( ) ( , ) pq pq c c
x x y y f x y dxdy

     . ( 2
)
For digital images, the above integral can be replaced by the following formula,

( ) ( ) ( , ) pq pq c c xy x x y y f x y      . ( 3
)
The Standardized central moment is

, pq  , which can be expressed mathematically as follows,

, , 0,0 2 , 2,3, 2 pq
pq pq where p q

         . (4)
In 1962, Hu [START_REF] Hu | Visual pattern recognition by moment invariant[END_REF] proposed the invariant moment theory of image recognition, and established a statistical feature extraction method for image recognition which has been widely used. Li [START_REF] Li | Reforming the theory of invariant moment for pattern recognition[END_REF] provided a calculation of 52 invariant moments, and seven common invariant moments are presented. From the above equation ( 1), (2) we could get seven u H moment invariants, which can be expressed as (3

1 20 02     . ( 5
) 22 2 20 02 11 ( ) 4        . ( 6
) 22 3 30 12 21 03 ( 3 ) (3 )          . ( 7
( 3 )( ) ( ) 3( ) (3 )( ) 3( ) ( )                                  . ( 9
                      . ( 10
)
)( ) ( ) 3( ) ( )( ) 3( ) ( )                                  . ( 11
)
Although the seven invariant moments can well describe shape character, but when the image database is large, only 7 scalar features are not enough for feature extraction. Because the information carried by their own is very limited, they usually need to be combined with other conventional and global shape features in order to well describe the actual shape properties of the object.

Definition of New Shape Feature

Seven u H invariant moments were selected as shape features in this paper, while area, perimeter, density and other parameters were selected as initial features. And two parameters were defined as shape features in intelligent identification system, which are as follows, 1) Number of spots (NS), which mainly used to calculate the number of spots in a disease leaf. 2) Ratio of area and number of spots (RANS). This parameter is an area of measurement of single spot, which mainly used to distinguish between large and small lesion spots.

Result and Discussion

In this study two categories of wheat rust and wheat leaf blight are as the main object, and the features are calculated according to the above definition formula. The effect drawing of shape feature extraction is shown in Figure 1.

Fig. 1. Interface of shape feature extraction

From experimental statistics we found that the above-defined parameters of moment invariants have a good translation, scale and rotation transformation invariance. While the defined three parameters of ratio of area and number of spots, density and number of spots can be identified as core features of identification of wheat rust and wheat leaf blight. When the value of ratio of area and number of spots is more than 1200, the wheat disease can be considered as the leaf blight disease. And when the parameter is less than 500, the disease can be considered as wheat rust disease. For the parameters of density and number of spots, they need to be used in combination with the parameter of ratio of area and number of spots. When the parameter of ratio of area and number of spots is used as identification basis, the recognition rate of wheat diseases is about 95.3%. When the three parameters of ratio of area and number of spots, density and number of spots are used as identification basis, the recognition rate of wheat diseases is about 97.1%. When invariant moments are taken into consider, and support vector machine algorithm is used as the model of disease identification, the recognition rate of wheat diseases is up to 99.4%. Using the above-defined shape features as the input of wheat disease identification system, we could obtain a higher recognition efficiency which can fully meet the needs of practical application.

Conclusion

Feature extraction process is a process to reduce the entropy of image information, and how to define and use features to express the image is the basis of image recognition [START_REF] Li | Reforming the theory of invariant moment for pattern recognition[END_REF]. In this study, combined with the current research of shape feature extraction, the invariant moment theory was introduced to shape feature extraction and seven Hu invariant moment parameters were defined as shape features. Meanwhile new parameters were defined for shape feature extraction research on wheat leaf disease image. The results showed that the system recognition rate is relatively high, and can meet the practical application requirements.

The shape features in this study is relatively easy to extract, and could reflect the morphological characteristics of two kinds of wheat disease images. But how to extract more features that can directly reflect the characteristics of different wheat disease images is the direction of future efforts.

  )

					(   	22 ) ( )    	.	(8)
			4		30	12		21	03
									22
	5	30	12		30	12	30	12		21	03
								22
	21	03		21	03		30	12	21	03

  )

					22
	(		) (		) (		)	4 (		)(		)
	6	20	02	30	12	21	03	11	30	12	21	03

This work was financially supported by the National High Tech. R&D Program of China (2007AA10Z237), and Doctor Foundation Project of Zhengzhou University of Light Industry. Author: Diao Zhihua(1982-), male, lecturer, doctor of engineering, interesting researches are image processing and multi-class classification algorithm.

Shape Feature ExtractionShape feature is an important feature to describe the image. Firstly, shape often links to the target, which has a certain semantic meaning. So shape feature can be seen higher level feature than color and texture feature. In addition, the shape expression is much more complex than the expression of the color or texture in nature. Secondly, the description of the target shape is a very complex issue. In fact, we still have not found the exact mathematical definition of shapes, including geometric, statistical, or morphological definition, so that it can be consistent with people's feelings. People's feeling of shape is a general result between retina feelings and knowledge about the real world of people. In order to accurately match the shape, the problem of translation, scaling, rotation invariance needs to be solved. As the target shape of images obtained from different vision may be involved in a large difference, translation, scale and rotation transformation invariance problems need to be solved for accurate shape matching.

Acknowledgments. This work was financially supported by the National High Tech. R&D Program of China (2007AA10Z237), and Doctor Foundation Project of Zhengzhou University of Light Industry.