Index: /issm/trunk-jpl/src/c/Container/Observations.cpp
===================================================================
--- /issm/trunk-jpl/src/c/Container/Observations.cpp	(revision 12297)
+++ /issm/trunk-jpl/src/c/Container/Observations.cpp	(revision 12298)
@@ -36,8 +36,9 @@
 
 	/*Intermediaries*/
-	int          i,maxdepth,level,counter;
+	int          i,j,maxdepth,level,counter,index;
 	int          xi,yi;
 	double       xmin,xmax,ymin,ymax;
-	double       offset,minlength;
+	double       offset,minlength,minspacing,mintrimming,maxtrimming;
+	int         *indices     = NULL;
 	Observation *observation = NULL;
 
@@ -52,4 +53,10 @@
 	offset=0.05*(ymax-ymin); ymin-=offset; ymax+=offset;
 
+	/*Get trimming limits*/
+	options->Get(&mintrimming,"mintrimming",-1.e+21);
+	options->Get(&maxtrimming,"maxtrimming",+1.e+21);
+	options->Get(&minspacing,"minspacing",0.01);
+	if(minspacing<=0) _error_("minspacing must > 0");
+
 	/*Get Minimum box size*/
 	if(options->GetOption("boxlength")){
@@ -70,4 +77,16 @@
 	counter = 0;
 	for(i=0;i<n;i++){
+
+		/*First check limits*/
+		if(observations_list[i]>maxtrimming) continue;
+		if(observations_list[i]<mintrimming) continue;
+
+		/*First check that this observation is not too close from another one*/
+		this->quadtree->ClosestObs(&index,x[i],y[i]);
+		if(index>=0){
+			observation=(Observation*)this->GetObjectByOffset(index);
+			if(pow(observation->x-x[i],2)+pow(observation->y-y[i],2) < minspacing) continue;
+		}
+
 		this->quadtree->IntergerCoordinates(&xi,&yi,x[i],y[i]);
 		this->quadtree->QuadtreeDepth2(&level,xi,yi);
@@ -83,4 +102,6 @@
 	}
 	printf("done\n");
+	printf("Initial number of observations: %i\n",n);
+	printf("  Final number of observations: %i\n",this->quadtree->NbObs);
 }
 /*}}}*/
@@ -98,7 +119,8 @@
 	/*Output and Intermediaries*/
 	bool         stop;
-	int          nobs,i,j,k,n,counter;
+	int          nobs,tempnobs,i,j,k,n,counter;
 	double       h2,radius2;
 	int         *indices      = NULL;
+	int         *tempindices  = NULL;
 	double      *dists        = NULL;
 	double      *x            = NULL;
@@ -112,15 +134,17 @@
 
 	/*Find all observations that are in radius*/
-	indices = (int*)xmalloc(this->Size()*sizeof(int));
-	dists   = (double*)xmalloc(this->Size()*sizeof(double));
-	nobs    = 0;
-
-	for (i=0;i<this->Size();i++){
-		observation=(Observation*)this->GetObjectByOffset(i);
+	this->quadtree->RangeSearch(&tempindices,&tempnobs,x_interp,y_interp,radius);
+	if(tempnobs){
+		indices = (int*)xmalloc(tempnobs*sizeof(int));
+		dists   = (double*)xmalloc(tempnobs*sizeof(double));
+	}
+	nobs = 0;
+	for (i=0;i<tempnobs;i++){
+		observation=(Observation*)this->GetObjectByOffset(tempindices[i]);
 		h2 = (observation->x-x_interp)*(observation->x-x_interp) + (observation->y-y_interp)*(observation->y-y_interp);
 
 		if(nobs==maxdata && h2>radius2) continue;
 		if(nobs<=maxdata){
-			indices[nobs]   = i;
+			indices[nobs]   = tempindices[i];
 			dists[nobs]     = h2;
 			nobs++;
@@ -141,5 +165,5 @@
 				}
 				dists[k]   = h2;
-				indices[k] = i;
+				indices[k] = tempindices[i];
 				stop = true;
 				break;
@@ -149,4 +173,5 @@
 	}  
 	xfree((void**)&dists);
+	xfree((void**)&tempindices);
 
 	if(nobs){
Index: /issm/trunk-jpl/src/c/objects/Kriging/Quadtree.cpp
===================================================================
--- /issm/trunk-jpl/src/c/objects/Kriging/Quadtree.cpp	(revision 12297)
+++ /issm/trunk-jpl/src/c/objects/Kriging/Quadtree.cpp	(revision 12298)
@@ -250,4 +250,47 @@
 	}
 }/*}}}*/
+/*FUNCTION Quadtree::ClosestObs{{{1*/
+void Quadtree::ClosestObs(int *pindex,double x,double y){
+
+	QuadtreeBox **pbox = NULL;
+	QuadtreeBox  *box  = NULL;
+	int           xi,yi;
+	int           level,levelbin;
+	int           index = -1;
+	double        length,length2;
+
+	/*Get integer coodinates*/
+	this->IntergerCoordinates(&xi,&yi,x,y);
+
+	/*Initialize levels*/
+	level    = 0;
+	levelbin = (1L<<this->MaxDepth);// = 2^30
+
+	/*Get inital box (the largest)*/
+	pbox=&root;
+
+	/*Find the smallest box where this point is located*/
+	while((box=*pbox) && (box->nbitems<0)){ 
+
+		levelbin>>=1; level+=1; 
+
+		pbox = &box->box[IJ(xi,yi,levelbin)];
+	}
+
+	/*Add obervation in this box (should be full)*/
+	if(box && box->nbitems>0){
+		index  = box->obs[0]->index;
+		length = pow(box->obs[0]->x - x,2.) + pow(box->obs[0]->y - y,2.);
+		for(int i=1;i<box->nbitems;i++){
+			length2 = pow(box->obs[i]->x - x,2.) + pow(box->obs[i]->y - y,2.);
+			if(length2<length){
+				index  = box->obs[i]->index;
+				length = length2;
+			}
+		}
+	}
+
+	*pindex=index;
+}/*}}}*/
 /*FUNCTION Quadtree::Echo{{{1*/
 void  Quadtree::Echo(void){
@@ -453,8 +496,8 @@
 
 	/*Allocate indices (maximum by default*/
-	indices = (int*)xmalloc(this->NbObs*sizeof(int));
+	if(this->NbObs) indices = (int*)xmalloc(this->NbObs*sizeof(int));
 	nobs = 0;
 
-	this->root->RangeSearch(indices,&nobs,x,y,range);
+	if(this->root) this->root->RangeSearch(indices,&nobs,x,y,range);
 
 	/*Clean-up and return*/
Index: /issm/trunk-jpl/src/c/objects/Kriging/Quadtree.h
===================================================================
--- /issm/trunk-jpl/src/c/objects/Kriging/Quadtree.h	(revision 12297)
+++ /issm/trunk-jpl/src/c/objects/Kriging/Quadtree.h	(revision 12298)
@@ -34,7 +34,7 @@
 
 				/*Methods*/
-				int IsWithinRange(double x,double y,double range);
-				void RangeSearch(int* indices,int *pnobs,double x,double y,double range);
-				void WriteObservations(int* indices,int *pnobs);
+				int          IsWithinRange(double  x,double y,double range);
+				void         RangeSearch(int *indices,int *pnobs,double x,double y,double range);
+				void         WriteObservations(int *indices,int *pnobs);
 
 		};
@@ -54,4 +54,5 @@
 		void         Add(Observation *observation);
 		void         AddAndAverage(double x,double y,double value);
+		void         ClosestObs(int *pindex,double x,double y);
 		void         DeepEcho(void);
 		void         Echo(void);
Index: /issm/trunk-jpl/src/modules/Kriging/Kriging.cpp
===================================================================
--- /issm/trunk-jpl/src/modules/Kriging/Kriging.cpp	(revision 12297)
+++ /issm/trunk-jpl/src/modules/Kriging/Kriging.cpp	(revision 12298)
@@ -55,5 +55,19 @@
 void KrigingUsage(void){
 	_printf_(true,"\n");
-	_printf_(true,"   usage: predictions=%s(x,y,observations,x_interp,y_interp);\n",__FUNCT__);
+	_printf_(true,"   usage: predictions=%s(x,y,observations,x_interp,y_interp,'options');\n",__FUNCT__);
+	_printf_(true,"   available options:\n");
+	_printf_(true,"      -'model': Available variogram models 'gaussian' (default),'spherical','power','exponential'\n");
+	_printf_(true,"         -'nugget': nugget effect (default 0.2)\n");
+	_printf_(true,"         -'range':  for gaussian, spherical and exponential models (default sqrt(3))\n");
+	_printf_(true,"         -'sill':   for gaussian, spherical and exponential models (default 1)\n");
+	_printf_(true,"         -'slope':  for power model (default 1)\n");
+	_printf_(true,"         -'power':  for power model (default 1)\n");
+	_printf_(true,"      -'searchradius': search radius for each prediction (default is observations span)\n");
+	_printf_(true,"      -'boxlength':    minimum length of quadtree boxes (useful to decrease the number of observations)\n");
+	_printf_(true,"      -'maxdata':      minimum number of observations for a prediction (default is 50)\n");
+	_printf_(true,"      -'mindata':      maximum number of observations for a prediction (default is 1)\n");
+	_printf_(true,"      -'maxtrimming':  maximum trimming value (default is -1.e+21)\n");
+	_printf_(true,"      -'mintrimming':  minimum trimming value (default is +1.e+21)\n");
+	_printf_(true,"      -'minspacing':   minimum distance between observation (default is 0.01)\n");
 	_printf_(true,"\n");
 }
