Precise GIS functions added.

12 files changed, 4898 insertions, 147 deletions

diff --git a/sql/CMakeLists.txt b/sql/CMakeLists.txt
index 8d2aeca17db..96c3643ccbc 100644
--- a/sql/CMakeLists.txt
+++ b/sql/CMakeLists.txt
@@ -78,6 +78,7 @@ SET (SQL_SOURCE
                rpl_rli.cc rpl_mi.cc sql_servers.cc
                sql_connect.cc scheduler.cc 
                sql_profile.cc event_parse_data.cc opt_table_elimination.cc
+               gcalc_slicescan.cc gcalc_tools.cc
                multi_range_read.cc
                opt_subselect.cc
                opt_index_cond_pushdown.cc
diff --git a/sql/Makefile.am b/sql/Makefile.am
index cde8962b8d0..54a2fafb2c7 100644
--- a/sql/Makefile.am
+++ b/sql/Makefile.am
@@ -56,6 +56,7 @@ noinst_HEADERS =	item.h item_func.h item_sum.h item_cmpfunc.h \
 			sql_map.h sql_string.h unireg.h \
 			sql_error.h field.h handler.h mysqld_suffix.h \
 			sql_profile.h \
+			gcalc_slicescan.h gcalc_tools.h \
 			ha_ndbcluster.h ha_ndbcluster_cond.h \
 			ha_ndbcluster_binlog.h ha_ndbcluster_tables.h \
 			ha_partition.h rpl_constants.h \
@@ -101,6 +102,7 @@ mysqld_SOURCES =	sql_lex.cc sql_handler.cc sql_partition.cc \
                         sql_join_cache.cc \
 			sql_base.cc table.cc sql_select.cc sql_insert.cc \
 			sql_profile.cc \
+			gcalc_slicescan.cc gcalc_tools.cc \
 			sql_prepare.cc sql_error.cc sql_locale.cc \
 			sql_update.cc sql_delete.cc uniques.cc sql_do.cc \
 			procedure.cc sql_test.cc \
diff --git a/sql/gcalc_slicescan.cc b/sql/gcalc_slicescan.cc
new file mode 100644
index 00000000000..4a77f309f46
--- /dev/null
+++ b/sql/gcalc_slicescan.cc
@@ -0,0 +1,762 @@
+/* Copyright (c) 2000, 2010 Oracle and/or its affiliates. All rights reserved.
+
+   This program is free software; you can redistribute it and/or modify
+   it under the terms of the GNU General Public License as published by
+   the Free Software Foundation; version 2 of the License.
+
+   This program is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+   GNU General Public License for more details.
+
+   You should have received a copy of the GNU General Public License
+   along with this program; if not, write to the Free Software
+   Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA */
+
+
+#include "mysql_priv.h"
+
+#ifdef HAVE_SPATIAL
+
+#include "gcalc_slicescan.h"
+
+
+#define PH_DATA_OFFSET 8
+#define coord_to_float(d) ((double) d)
+
+typedef int (*sc_compare_func)(const void*, const void*);
+
+#define LS_LIST_ITEM Gcalc_dyn_list::Item
+#define LS_COMPARE_FUNC_DECL sc_compare_func compare,
+#define LS_COMPARE_FUNC_CALL(list_el1, list_el2) (*compare)(list_el1, list_el2)
+#define LS_NEXT(A) (A)->next
+#define LS_SET_NEXT(A,val) (A)->next= val
+#define LS_P_NEXT(A) &(A)->next
+#define LS_NAME sort_list
+#define LS_SCOPE static
+#define LS_STRUCT_NAME sort_list_stack_struct
+#include "plistsort.c"
+
+
+Gcalc_dyn_list::Gcalc_dyn_list(size_t blk_size, size_t sizeof_item):
+  m_blk_size(blk_size - ALLOC_ROOT_MIN_BLOCK_SIZE),
+  m_sizeof_item(ALIGN_SIZE(sizeof_item)),
+  m_points_per_blk((m_blk_size - PH_DATA_OFFSET) / m_sizeof_item),
+  m_blk_hook(&m_first_blk),
+  m_free(NULL),
+  m_keep(NULL)
+{}
+
+
+void Gcalc_dyn_list::format_blk(void* block)
+{
+  Item *pi_end, *cur_pi, *first_pi;
+  DBUG_ASSERT(m_free == NULL);
+  first_pi= cur_pi= (Item *)(((char *)block) + PH_DATA_OFFSET);
+  pi_end= ptr_add(first_pi, m_points_per_blk - 1);
+  do {
+    cur_pi= cur_pi->next= ptr_add(cur_pi, 1);
+  } while (cur_pi<pi_end);
+  cur_pi->next= m_free;
+  m_free= first_pi;
+}
+
+
+Gcalc_dyn_list::Item *Gcalc_dyn_list::alloc_new_blk()
+{
+  void *new_block= my_malloc(m_blk_size, MYF(MY_WME));
+  if (!new_block)
+    return NULL;
+  *m_blk_hook= new_block;
+  m_blk_hook= (void**)new_block;
+  format_blk(new_block);
+  return new_item();
+}
+
+
+static void free_blk_list(void *list)
+{
+  void *next_blk;
+  while (list)
+  {
+    next_blk= *((void **)list);
+    my_free(list, MYF(0));
+    list= next_blk;
+  }
+}
+
+
+void Gcalc_dyn_list::cleanup()
+{
+  *m_blk_hook= NULL;
+  free_blk_list(m_first_blk);
+  m_first_blk= NULL;
+  m_blk_hook= &m_first_blk;
+  m_free= NULL;
+}
+
+
+Gcalc_dyn_list::~Gcalc_dyn_list()
+{
+  cleanup();
+}
+
+
+void Gcalc_dyn_list::reset()
+{
+  *m_blk_hook= NULL;
+  if (m_first_blk)
+  {
+    free_blk_list(*((void **)m_first_blk));
+    m_blk_hook= (void**)m_first_blk;
+    m_free= NULL;
+    format_blk(m_first_blk);
+  }
+}
+
+
+static inline void trim_node(Gcalc_heap::Info *node, Gcalc_heap::Info *prev_node)
+{
+  if (!node)
+    return;
+  DBUG_ASSERT((node->left == prev_node) || (node->right == prev_node));
+  if (node->left == prev_node)
+    node->left= node->right;
+  node->right= NULL;
+}
+
+
+static double find_first_different(const Gcalc_heap::Info *p)
+{
+  if (p->left && (p->left->y != p->y))
+    return p->left->y;
+  if (p->right && (p->right->y != p->y))
+    return p->right->y;
+  if (p->left && p->left->left && (p->left->left->y != p->y))
+    return p->left->left->y;
+  if (p->right && p->right->right && (p->right->right->y != p->y))
+    return p->right->right->y;
+
+  return p->y;
+}
+
+
+static int compare_point_info(const void *e0, const void *e1)
+{
+  const Gcalc_heap::Info *i0= (const Gcalc_heap::Info *)e0;
+  const Gcalc_heap::Info *i1= (const Gcalc_heap::Info *)e1;
+  if (i0->y != i1->y)
+    return i0->y > i1->y;
+  return find_first_different(i0) > find_first_different(i1);
+}
+
+
+void Gcalc_heap::prepare_operation()
+{
+  DBUG_ASSERT(m_hook);
+  *m_hook= NULL;
+  m_first= sort_list(compare_point_info, m_first, m_n_points);
+  m_hook= NULL; /* just to check it's not called twice */
+
+  /* TODO - move this to the 'normal_scan' loop */
+  for (Info *cur= get_first(); cur; cur= cur->get_next())
+  {
+    trim_node(cur->left, cur);
+    trim_node(cur->right, cur);
+  }
+}
+
+
+void Gcalc_heap::reset()
+{
+  if (!m_hook)
+  {
+    m_hook= &m_first;
+    for (; *m_hook; m_hook= &(*m_hook)->next)
+    {}
+  }
+
+  *m_hook= m_free;
+  m_free= m_first;
+  m_hook= &m_first;
+  m_n_points= 0;
+}
+
+int Gcalc_shape_transporter::int_single_point(gcalc_shape_info Info,
+                                              double x, double y)
+{
+  Gcalc_heap::Info *point= m_heap->new_point_info(x, y, Info);
+  if (!point)
+    return 1;
+  point->left= point->right= 0;
+  return 0;
+}
+
+
+int Gcalc_shape_transporter::int_add_point(gcalc_shape_info Info,
+                                           double x, double y)
+{
+  Gcalc_heap::Info *point;
+  if (!(point= m_heap->new_point_info(x, y, Info)))
+    return 1;
+  if (m_first)
+  {
+    m_prev->left= point;
+    point->right= m_prev;
+  }
+  else
+    m_first= point;
+  m_prev= point;
+  return 0;
+}
+
+
+void Gcalc_shape_transporter::int_complete()
+{
+  DBUG_ASSERT(m_shape_started == 1 || m_shape_started == 3);
+
+  if (!m_first)
+    return;
+
+  /* simple point */
+  if (m_first == m_prev)
+  {
+    m_first->right= m_first->left= NULL;
+    return;
+  }
+
+  /* line */
+  if (m_shape_started == 1)
+  {
+    m_first->right= NULL;
+    m_prev->left= m_prev->right;
+    m_prev->right= NULL;
+    return;
+  }
+
+  /* polygon */
+  m_first->right= m_prev;
+  m_prev->left= m_first;
+}
+
+
+inline int GET_DX_DY(double *dxdy,
+                     const Gcalc_heap::Info *p0, const Gcalc_heap::Info *p1)
+{
+  double dy= p1->y - p0->y;
+  *dxdy= p1->x - p0->x;
+  return (dy == 0.0) ||
+         (*dxdy/= dy)>DBL_MAX ||
+         (*dxdy)<-DBL_MAX;
+}
+
+Gcalc_scan_iterator::Gcalc_scan_iterator(size_t blk_size) :
+  Gcalc_dyn_list(blk_size,
+	         (sizeof(point) > sizeof(intersection)) ?
+	          sizeof(point) : sizeof(intersection)),
+  m_slice0(NULL), m_slice1(NULL)
+{}
+		  
+Gcalc_scan_iterator::point
+  *Gcalc_scan_iterator::new_slice(Gcalc_scan_iterator::point *example)
+{
+  point *result= NULL;
+  Gcalc_dyn_list::Item **result_hook= (Gcalc_dyn_list::Item **)&result;
+  while (example)
+  {
+    *result_hook= new_slice_point();
+    result_hook= &(*result_hook)->next;
+    example= example->get_next();
+  }
+  *result_hook= NULL;
+  return result;
+}
+
+
+void Gcalc_scan_iterator::init(Gcalc_heap *points)
+{
+  DBUG_ASSERT(points->ready());
+  DBUG_ASSERT(!m_slice0 && !m_slice1);
+
+  if (!(m_cur_pi= points->get_first()))
+    return;
+  m_cur_thread= 0;
+  m_sav_slice= NULL;
+  m_intersections= NULL;
+  m_cur_intersection= NULL;
+  m_y1= m_cur_pi->y;
+  m_next_is_top_point= true;
+  m_bottom_points_count= 0;
+}
+
+void Gcalc_scan_iterator::reset()
+{
+  if (m_slice0)
+    free_list(m_slice0);
+  if (m_slice1)
+    free_list(m_slice1);
+  m_slice0= m_slice1= NULL;
+  Gcalc_dyn_list::reset();
+}
+
+static bool slice_first_equal_x(const Gcalc_scan_iterator::point *p0,
+				const Gcalc_scan_iterator::point *p1)
+{
+  if (p0->horiz_dir == p1->horiz_dir)
+    return p0->dx_dy <= p1->dx_dy;
+  if (p0->horiz_dir)
+    return p0->dx_dy < 0;
+  return p1->dx_dy > 0;  /* p1->horiz_dir case */
+}
+
+
+static inline bool slice_first(const Gcalc_scan_iterator::point *p0,
+			       const Gcalc_scan_iterator::point *p1)
+{
+  if (p0->x != p1->x)
+    return p0->x < p1->x;
+  return slice_first_equal_x(p0, p1);
+}
+
+
+int Gcalc_scan_iterator::insert_top_point()
+{
+  point *sp= m_slice1;
+  Gcalc_dyn_list::Item **prev_hook= (Gcalc_dyn_list::Item **)&m_slice1;
+  point *sp1;
+  point *sp0= new_slice_point();
+
+  if (!sp0)
+    return 1;
+  sp0->pi= m_cur_pi;
+  sp0->next_pi= m_cur_pi->left;
+  sp0->thread= m_cur_thread++;
+  sp0->x= coord_to_float(m_cur_pi->x);
+  if (m_cur_pi->left)
+  {
+    sp0->horiz_dir= GET_DX_DY(&sp0->dx_dy, m_cur_pi, m_cur_pi->left);
+    m_event1= scev_thread;
+
+    /*Now just to increase the size of m_slice0 to be same*/
+    if (!(sp1= new_slice_point()))
+      return 1;
+    sp1->next= m_slice0;
+    m_slice0= sp1;
+  }
+  else
+  {
+    m_event1= scev_single_point;
+    sp0->horiz_dir= 0;
+    sp0->dx_dy= 0.0;
+  }
+
+  /* First we need to find the place to insert.
+     Binary search could probably make things faster here,
+     but structures used aren't suitable, and the
+     scan is usually not really long */
+  for (; sp && slice_first(sp, sp0);
+       prev_hook= &sp->next, sp=sp->get_next())
+  {}
+
+  if (m_cur_pi->right)
+  {
+    m_event1= scev_two_threads;
+    /*We have two threads so should decide which one will be first*/
+    sp1= new_slice_point();
+    if (!sp1)
+      return 1;
+    sp1->pi= m_cur_pi;
+    sp1->next_pi= m_cur_pi->right;
+    sp1->thread= m_cur_thread++;
+    sp1->x= sp0->x;
+    sp1->horiz_dir= GET_DX_DY(&sp1->dx_dy, m_cur_pi, m_cur_pi->right);
+    if (slice_first_equal_x(sp1, sp0))
+    {
+      point *tmp= sp0;
+      sp0= sp1;
+      sp1= tmp;
+    }
+    sp1->next= sp;
+    sp0->next= sp1;
+    
+    /*Now just to increase the size of m_slice0 to be same*/
+    if (!(sp1= new_slice_point()))
+      return 1;
+    sp1->next= m_slice0;
+    m_slice0= sp1;
+  }
+  else
+    sp0->next= sp;
+
+  *prev_hook= sp0;
+  m_event_position1= sp0;
+
+  return 0;
+}
+
+enum
+{
+  intersection_normal= 1,
+  intersection_forced= 2
+};
+
+
+static int intersection_found(const Gcalc_scan_iterator::point *sp0,
+			      const Gcalc_scan_iterator::point *sp1,
+			      unsigned int bottom_points_count)
+{
+  if (sp1->x < sp0->x)
+    return intersection_normal;
+  if (sp1->is_bottom() && !sp0->is_bottom() &&
+      (bottom_points_count > 1))
+      return intersection_forced;
+  return 0;
+}
+
+
+int Gcalc_scan_iterator::normal_scan()
+{
+  if (m_next_is_top_point)
+    if (insert_top_point())
+      return 1;
+
+  point *tmp= m_slice0;
+  m_slice0= m_slice1;
+  m_slice1= tmp;
+  m_event0= m_event1;
+  m_event_position0= m_event_position1;
+  m_y0= m_y1;
+  
+  if (!(m_cur_pi= m_cur_pi->get_next()))
+  {
+    free_list(m_slice1);
+    m_slice1= NULL;
+    return 0;
+  }
+  
+  Gcalc_heap::Info *cur_pi= m_cur_pi;
+  m_y1= coord_to_float(cur_pi->y);
+  m_h= m_y1 - m_y0;
+
+  point *sp0= m_slice0;
+  point *sp1= m_slice1;
+  point *prev_sp1= NULL;
+
+  m_bottom_points_count= 0;
+  m_next_is_top_point= true;
+  bool intersections_found= false;
+
+  for (; sp0; sp0= sp0->get_next())
+  {
+    if (sp0->next_pi == cur_pi) /* End of the segment */
+    {
+      sp1->x= coord_to_float(cur_pi->x);
+      sp1->pi= cur_pi;
+      sp1->thread= sp0->thread;
+      sp1->next_pi= cur_pi->left;
+      if (cur_pi->left)
+	sp1->horiz_dir= GET_DX_DY(&sp1->dx_dy, m_cur_pi, m_cur_pi->left);
+
+      m_next_is_top_point= false;
+      
+      if (sp1->is_bottom())
+      {
+	++m_bottom_points_count;
+	if (m_bottom_points_count == 1)
+	{
+	  m_event1= scev_end;
+	  m_event_position1= sp1;
+	}
+	else
+	  m_event1= scev_two_ends;
+      }
+      else
+      {
+	m_event1= scev_point;
+	m_event_position1= sp1;
+      }
+    }
+    else if (!sp0->is_bottom())
+    {
+      /* Cut current string with the height of the new point*/
+      sp1->copy_core(sp0);
+      sp1->x= sp1->horiz_dir ? sp0->x :
+	(coord_to_float(sp1->pi->x) +
+	 (m_y1-coord_to_float(sp1->pi->y)) * sp1->dx_dy);
+    }
+    else  /* Skip the bottom point in slice0 */
+      continue;
+
+    intersections_found= intersections_found ||
+      (prev_sp1 && intersection_found(prev_sp1, sp1, m_bottom_points_count));
+
+    prev_sp1= sp1;
+    sp1= sp1->get_next();
+  }
+
+  if (sp1)
+  {
+    if (prev_sp1)
+      prev_sp1->next= NULL;
+    else
+      m_slice1= NULL;
+    free_list(sp1);
+  }
+
+  if (intersections_found)
+    return handle_intersections();
+
+  return 0;
+}
+
+
+int Gcalc_scan_iterator::add_intersection(const point *a, const point *b,
+				   int isc_kind, Gcalc_dyn_list::Item ***p_hook)
+{
+  intersection *isc= new_intersection();
+
+  if (!isc)
+    return 1;
+  m_n_intersections++;
+  **p_hook= isc;
+  *p_hook= &isc->next;
+  isc->thread_a= a->thread;
+  isc->thread_b= b->thread;
+  if (isc_kind == intersection_forced)
+  {
+    isc->y= m_y1;
+    isc->x= a->x;
+    return 0;
+  }
+
+  /* intersection_normal */
+  const point *a0= a->precursor;
+  const point *b0= b->precursor;
+  if (!a0->horiz_dir && !b0->horiz_dir)
+  {
+    double dk= a0->dx_dy - b0->dx_dy;
+    double dy= (b0->x - a0->x)/dk;
+    isc->y= m_y0 + dy;
+    isc->x= a0->x + dy*a0->dx_dy;
+    return 0;
+  }
+  isc->y= m_y1;
+  isc->x= a0->horiz_dir ? b->x : a->x;
+  return 0;
+}
+
+
+int Gcalc_scan_iterator::find_intersections()
+{
+  point *sp1= m_slice1;
+  Gcalc_dyn_list::Item **hook;
+
+  m_n_intersections= 0;
+  {
+    /* Set links between slicepoints */
+    point *sp0= m_slice0;
+    for (; sp1; sp0= sp0->get_next(),sp1= sp1->get_next())
+    {
+      while (sp0->is_bottom())
+	sp0= sp0->get_next();
+      DBUG_ASSERT(sp0->thread == sp1->thread);
+      sp1->precursor= sp0;
+    }
+  }
+
+  hook= (Gcalc_dyn_list::Item **)&m_intersections;
+  bool intersections_found;
+
+  point *last_possible_isc= NULL;
+  do
+  {
+    sp1= m_slice1;
+    point **pprev_s1= &m_slice1;
+    intersections_found= false;
+    unsigned int bottom_points_count= sp1->is_bottom() ? 1:0;
+    sp1= m_slice1->get_next();
+    int isc_kind;
+    point *cur_possible_isc= NULL;
+    for (; sp1 != last_possible_isc;
+	 pprev_s1= (point **)(&(*pprev_s1)->next), sp1= sp1->get_next())
+    {
+      if (sp1->is_bottom())
+	++bottom_points_count;
+      if (!(isc_kind=intersection_found(*pprev_s1, sp1, bottom_points_count)))
+	continue;
+      point *prev_s1= *pprev_s1;
+      intersections_found= true;
+      if (add_intersection(prev_s1, sp1, isc_kind, &hook))
+	return 1;
+      *pprev_s1= sp1;
+      prev_s1->next= sp1->next;
+      sp1->next= prev_s1;
+      sp1= prev_s1;
+      cur_possible_isc= sp1;
+    }
+    last_possible_isc= cur_possible_isc;
+  } while (intersections_found);
+
+  *hook= NULL;
+  return 0;
+}
+
+
+static int compare_intersections(const void *e0, const void *e1)
+{
+  Gcalc_scan_iterator::intersection *i0= (Gcalc_scan_iterator::intersection *)e0;
+  Gcalc_scan_iterator::intersection *i1= (Gcalc_scan_iterator::intersection *)e1;
+  return i0->y > i1->y;
+}
+
+
+inline void Gcalc_scan_iterator::sort_intersections()
+{
+  m_intersections= (intersection *)sort_list(compare_intersections,
+                                             m_intersections,m_n_intersections);
+}
+
+
+int Gcalc_scan_iterator::handle_intersections()
+{
+  DBUG_ASSERT(m_slice1->next);
+
+  if (find_intersections())
+    return 1;
+  sort_intersections();
+
+  m_sav_slice= m_slice1;
+  m_sav_y= m_y1;
+  m_slice1= new_slice(m_sav_slice);
+  
+  m_cur_intersection= m_intersections;
+  m_pre_intersection_hook= NULL;
+  return intersection_scan();
+}
+
+
+void Gcalc_scan_iterator::pop_suitable_intersection()
+{
+  intersection *prev_i= m_cur_intersection;
+  intersection *cur_i= prev_i->get_next();
+  for (; cur_i; prev_i= cur_i, cur_i= cur_i->get_next())
+  {
+    point *prev_p= m_slice0;
+    point *sp= prev_p->get_next();
+    for (; sp; prev_p= sp, sp= sp->get_next())
+    {
+      if ((prev_p->thread == cur_i->thread_a) &&
+	  (sp->thread == cur_i->thread_b))
+      {
+	/* Move cur_t on the top of the list */
+	if (prev_i == m_cur_intersection)
+	{
+	  m_cur_intersection->next= cur_i->next;
+	  cur_i->next= m_cur_intersection;
+	  m_cur_intersection= cur_i;
+	}
+	else
+	{
+          Gcalc_dyn_list::Item *tmp= m_cur_intersection->next;
+	  m_cur_intersection->next= cur_i->next;
+	  prev_i->next= m_cur_intersection;
+	  m_cur_intersection= cur_i;
+	  cur_i->next= tmp;
+	}
+	return;
+      }
+    }
+  }
+  DBUG_ASSERT(0);
+}
+
+
+int Gcalc_scan_iterator::intersection_scan()
+{
+  if (m_pre_intersection_hook) /*Skip the first point*/
+  {
+    point *next= (*m_pre_intersection_hook)->get_next();
+    (*m_pre_intersection_hook)->next= next->next;
+    next->next= *m_pre_intersection_hook;
+    *m_pre_intersection_hook= next;
+    m_event0= scev_intersection;
+    m_event_position0= next;
+    point *tmp= m_slice1;
+    m_slice1= m_slice0;
+    m_slice0= tmp;
+    m_y0= m_y1;
+    m_cur_intersection= m_cur_intersection->get_next();
+    if (!m_cur_intersection)
+    {
+      m_h= m_sav_y - m_y1;
+      m_y1= m_sav_y;
+      free_list(m_slice1);
+      m_slice1= m_sav_slice;
+      free_list(m_intersections);
+      return 0;
+    }
+  }
+
+  m_y1= m_cur_intersection->y;
+  m_h= m_y1 - m_y0;
+
+  point *sp0;
+  point **psp1;
+
+redo_loop:
+  sp0= m_slice0;
+  psp1= &m_slice1;
+  for (; sp0; sp0= sp0->get_next())
+  {
+    point *sp1= *psp1;
+    if (sp0->thread == m_cur_intersection->thread_a)
+    {
+      point *next_s0= sp0;
+      /* Skip Bottom points */
+      do
+	next_s0= next_s0->get_next();
+      while(next_s0->is_bottom()); /* We always find nonbottom point here*/
+      /* If the next point's thread isn't the thread of intersection,
+	 we try to find suitable intersection */
+      if (next_s0->thread != m_cur_intersection->thread_b)
+      {
+	/* It's really rare case - sometimes happen when
+	   there's two intersections with the same Y
+	   Move suitable one to the beginning of the list
+	*/
+	pop_suitable_intersection();
+	goto redo_loop;
+      }
+      m_pre_intersection_hook= psp1;
+      sp1->copy_core(sp0);
+      sp1->x= m_cur_intersection->x;
+      sp0= next_s0;
+      sp1= sp1->get_next();
+      sp1->copy_core(sp0);
+      sp1->x= m_cur_intersection->x;
+      psp1= (point **)&sp1->next;
+      continue;
+    }
+    if (!sp0->is_bottom())
+    {
+      sp1->copy_core(sp0);
+      sp1->x= sp1->horiz_dir ? sp0->x :
+	(coord_to_float(sp1->pi->x) +
+	 (m_y1-coord_to_float(sp1->pi->y)) * sp1->dx_dy);
+    }
+    else
+      /* Skip bottom point */
+      continue;
+    psp1= (point **)&sp1->next;
+  }
+
+  if (*psp1)
+  {
+    free_list(*psp1);
+    *psp1= NULL;
+  }
+
+  return 0;
+}
+
+#endif /* HAVE_SPATIAL */
diff --git a/sql/gcalc_slicescan.h b/sql/gcalc_slicescan.h
new file mode 100644
index 00000000000..3ea910dc4fd
--- /dev/null
+++ b/sql/gcalc_slicescan.h
@@ -0,0 +1,433 @@
+/* Copyright (c) 2000, 2010 Oracle and/or its affiliates. All rights reserved.
+
+   This program is free software; you can redistribute it and/or modify
+   it under the terms of the GNU General Public License as published by
+   the Free Software Foundation; version 2 of the License.
+
+   This program is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+   GNU General Public License for more details.
+
+   You should have received a copy of the GNU General Public License
+   along with this program; if not, write to the Free Software
+   Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA */
+
+
+#ifndef GCALC_SLICESCAN_INCLUDED
+#define GCALC_SLICESCAN_INCLUDED
+
+
+/*
+  Gcalc_dyn_list class designed to manage long lists of same-size objects
+  with the possible efficiency.
+  It allocates fixed-size blocks of memory (blk_size specified at the time
+  of creation). When new object is added to the list, it occupies part of
+  this block until it's full. Then the new block is allocated.
+  Freed objects are chained to the m_free list, and if it's not empty, the
+  newly added object is taken from this list instead the block.
+*/
+
+class Gcalc_dyn_list
+{
+public:
+  class Item
+  {
+  public:
+    Item *next;
+  };
+
+  Gcalc_dyn_list(size_t blk_size, size_t sizeof_item);
+  ~Gcalc_dyn_list();
+  Item *new_item()
+  {
+    Item *result;
+    if (m_free)
+    {
+      result= m_free;
+      m_free= m_free->next;
+    }
+    else
+      result= alloc_new_blk();
+
+    return result;
+  }
+  inline void free_item(Item *item)
+  {
+    item->next= m_free;
+    m_free= item;
+  }
+  inline void free_list(Item *list, Item **hook)
+  {
+    *hook= m_free;
+    m_free= list;
+  }
+
+  void free_list(Item *list)
+  {
+    Item **hook= &list;
+    while (*hook)
+      hook= &(*hook)->next;
+    free_list(list, hook);
+  }
+
+  void reset();
+  void cleanup();
+
+protected:
+  size_t m_blk_size;
+  size_t m_sizeof_item;
+  unsigned int m_points_per_blk;
+  void *m_first_blk;
+  void **m_blk_hook;
+  Item *m_free;
+  Item *m_keep;
+
+  Item *alloc_new_blk();
+  void format_blk(void* block);
+  inline Item *ptr_add(Item *ptr, int n_items)
+  {
+    return (Item *)(((char*)ptr) + n_items * m_sizeof_item);
+  }
+};
+
+typedef uint gcalc_shape_info;
+
+/*
+  Gcalc_heap represents the 'dynamic list' of Info objects, that
+  contain information about vertexes of all the shapes that take
+  part in some spatial calculation. Can become quite long.
+  After filled, the list is usually sorted and then walked through
+  in the slicescan algorithm.
+  The Gcalc_heap and the algorithm can only operate with two
+  kinds of shapes - polygon and polyline. So all the spatial
+  objects should be represented as sets of these two.
+*/
+
+class Gcalc_heap : public Gcalc_dyn_list
+{
+public:
+  class Info : public Gcalc_dyn_list::Item
+  {
+  public:
+    gcalc_shape_info shape;
+    Info *left;
+    Info *right;
+    double x,y;
+
+    inline bool is_bottom() const { return !left; }
+    inline Info *get_next() { return (Info *)next; }
+    inline const Info *get_next() const { return (const Info *)next; }
+  };
+
+  Gcalc_heap(size_t blk_size=8192) :
+    Gcalc_dyn_list(blk_size, sizeof(Info)), m_hook(&m_first), m_n_points(0) {}
+  Info *new_point_info(double x, double y, gcalc_shape_info shape)
+  {
+    Info *result= (Info *)new_item();
+    if (!result)
+      return NULL;
+    *m_hook= result;
+    m_hook= &result->next;
+    m_n_points++;
+    result->x= x;
+    result->y= y;
+    result->shape= shape;
+    return result;
+  }
+  void prepare_operation();
+  inline bool ready() const { return m_hook == NULL; }
+  Info *get_first() { return (Info *)m_first; }
+  const Info *get_first() const { return (const Info *)m_first; }
+  Gcalc_dyn_list::Item **get_last_hook() { return m_hook; }
+  void reset();
+private:
+  Gcalc_dyn_list::Item *m_first;
+  Gcalc_dyn_list::Item **m_hook;
+  int m_n_points;
+};
+
+
+/*
+  the spatial object has to be represented as a set of
+  simple polygones and polylines to be sent to the slicescan.
+
+  Gcalc_shape_transporter class and his descendants are used to
+  simplify storing the information about the shape into necessary structures.
+  This base class only fills the Gcalc_heap with the information about
+  shapes and vertices.
+
+  Normally the Gcalc_shape_transporter family object is sent as a parameter
+  to the 'get_shapes' method of an 'spatial' object so it can pass
+  the spatial information about itself. The virtual methods are
+  treating this data in a way the caller needs.
+*/
+
+class Gcalc_shape_transporter
+{
+private:
+  Gcalc_heap::Info *m_first;
+  Gcalc_heap::Info *m_prev;
+  int m_shape_started;
+  void int_complete();
+protected:
+  Gcalc_heap *m_heap;
+  int int_single_point(gcalc_shape_info Info, double x, double y);
+  int int_add_point(gcalc_shape_info Info, double x, double y);
+  void int_start_line()
+  {
+    DBUG_ASSERT(!m_shape_started);
+    m_shape_started= 1;
+    m_first= m_prev= NULL;
+  }
+  void int_complete_line()
+  {
+    DBUG_ASSERT(m_shape_started== 1);
+    int_complete();
+    m_shape_started= 0;
+  }
+  void int_start_ring()
+  {
+    DBUG_ASSERT(m_shape_started== 2);
+    m_shape_started= 3;
+    m_first= m_prev= NULL;
+  }
+  void int_complete_ring()
+  {
+    DBUG_ASSERT(m_shape_started== 3);
+    int_complete();
+    m_shape_started= 2;
+  }
+  void int_start_poly()
+  {
+    DBUG_ASSERT(!m_shape_started);
+    m_shape_started= 2;
+  }
+  void int_complete_poly()
+  {
+    DBUG_ASSERT(m_shape_started== 2);
+    m_shape_started= 0;
+  }
+  bool line_started() { return m_shape_started == 1; };
+public:
+  Gcalc_shape_transporter(Gcalc_heap *heap) :
+    m_shape_started(0), m_heap(heap) {}
+
+  virtual int single_point(double x, double y)=0;
+  virtual int start_line()=0;
+  virtual int complete_line()=0;
+  virtual int start_poly()=0;
+  virtual int complete_poly()=0;
+  virtual int start_ring()=0;
+  virtual int complete_ring()=0;
+  virtual int add_point(double x, double y)=0;
+  virtual int start_collection(int n_objects) { return 0; }
+  int start_simple_poly()
+  {
+    return start_poly() || start_ring();
+  }
+  int complete_simple_poly()
+  {
+    return complete_ring() || complete_poly();
+  }
+  virtual ~Gcalc_shape_transporter() {}
+};
+
+
+enum Gcalc_scan_events
+{
+  scev_point= 1,         /* Just a new point in thread */
+  scev_thread= 2,        /* Start of the new thread */
+  scev_two_threads= 4,   /* A couple of new threads started */
+  scev_intersection= 8,  /* Intersection happened */
+  scev_end= 16,          /* Single thread finished */
+  scev_two_ends= 32,     /* A couple of threads finished */
+  scev_single_point= 64  /* Got single point */
+};
+
+typedef int sc_thread_id;
+
+/* 
+   Gcalc_scan_iterator incapsulates the slisescan algorithm.
+   It takes filled Gcalc_heap as an datasource. Then can be
+   iterated trought the vertexes and intersection points with
+   the step() method. After the 'step()' one usually observes
+   the current 'slice' to do the necessary calculations, like
+   looking for intersections, calculating the area, whatever.
+*/
+
+class Gcalc_scan_iterator : public Gcalc_dyn_list
+{
+public:
+  class point : public Gcalc_dyn_list::Item
+  {
+  public:
+    double x;
+    double dx_dy;
+    int horiz_dir;
+    Gcalc_heap::Info *pi;
+    Gcalc_heap::Info *next_pi;
+    sc_thread_id thread;
+    const point *precursor; /* used as a temporary field */
+
+    inline const point *c_get_next() const
+      { return (const point *)next; }
+    inline bool is_bottom() const { return pi->is_bottom(); }
+    gcalc_shape_info get_shape() const { return pi->shape; }
+    inline point *get_next() { return (point *)next; }
+    inline const point *get_next() const { return (const point *)next; }
+
+    /* copies all but 'next' 'x' and 'precursor' */
+    void copy_core(const point *from)
+    {
+      dx_dy= from->dx_dy;
+      horiz_dir= from->horiz_dir;
+      pi= from->pi;
+      next_pi= from->next_pi;
+      thread= from->thread;
+    }
+#ifndef DBUG_OFF
+    void dbug_print();
+#endif /*DBUG_OFF*/
+  };
+
+  class intersection : public Gcalc_dyn_list::Item
+  {
+  public:
+    sc_thread_id thread_a;
+    sc_thread_id thread_b;
+    double x;
+    double y;
+    inline intersection *get_next() { return (intersection *)next; }
+  };
+
+public:
+  Gcalc_scan_iterator(size_t blk_size= 8192);
+
+  void init(Gcalc_heap *points); /* Iterator can be reused */
+  void reset();
+  int step()
+  {
+    DBUG_ASSERT(more_points());
+    return m_cur_intersection ? intersection_scan() : normal_scan();
+  }
+
+  inline Gcalc_heap::Info *more_points() { return m_cur_pi; }
+  inline bool more_trapezoids()
+    { return m_cur_pi && m_cur_pi->next; }
+
+  inline Gcalc_scan_events get_event() const { return m_event0; }
+  inline const point *get_event_position() const
+    { return m_event_position0; }
+  inline const point *get_b_slice() const { return m_slice0; }
+  inline const point *get_t_slice() const { return m_slice1; }
+  inline double get_h() const { return m_h; }
+  inline double get_y() const { return m_y0; }
+
+private:
+  Gcalc_heap::Info *m_cur_pi;
+  point *m_slice0;
+  point *m_slice1;
+  point *m_sav_slice;
+  intersection *m_intersections;
+  int m_n_intersections;
+  intersection *m_cur_intersection;
+  point **m_pre_intersection_hook;
+  double m_h;
+  double m_y0;
+  double m_y1;
+  double m_sav_y;
+  bool m_next_is_top_point;
+  unsigned int m_bottom_points_count;
+  sc_thread_id m_cur_thread;
+  Gcalc_scan_events m_event0, m_event1;
+  point *m_event_position0;
+  point *m_event_position1;
+
+  int normal_scan();
+  int intersection_scan();
+  void sort_intersections();
+  int handle_intersections();
+  int insert_top_point();
+  int add_intersection(const point *a, const point *b,
+		       int isc_kind, Gcalc_dyn_list::Item ***p_hook);
+  int find_intersections();
+  void pop_suitable_intersection();
+
+  intersection *new_intersection()
+  {
+    return (intersection *)new_item();
+  }
+  point *new_slice_point()
+  {
+    return (point *)new_item();
+  }
+  point *new_slice(point *example);
+};
+
+
+/* 
+   Gcalc_trapezoid_iterator simplifies the calculations on
+   the current slice of the Gcalc_scan_iterator.
+   One can walk through the trapezoids formed between
+   previous and current slices.
+*/
+
+class Gcalc_trapezoid_iterator
+{
+protected:
+  const Gcalc_scan_iterator::point *sp0;
+  const Gcalc_scan_iterator::point *sp1;
+public:
+  Gcalc_trapezoid_iterator(const Gcalc_scan_iterator *scan_i) :
+    sp0(scan_i->get_b_slice()),
+    sp1(scan_i->get_t_slice())
+    {}
+
+  inline bool more() const { return sp1 && sp1->next; }
+
+  const Gcalc_scan_iterator::point *lt() const { return sp1; }
+  const Gcalc_scan_iterator::point *lb() const { return sp0; }
+  const Gcalc_scan_iterator::point *rb() const
+  {
+    const Gcalc_scan_iterator::point *result= sp0;
+    while ((result= result->c_get_next())->is_bottom())
+    {}
+    return result;
+  }
+  const Gcalc_scan_iterator::point *rt() const
+    { return sp1->c_get_next(); }
+
+  void operator++()
+  {
+    sp0= rb();
+    sp1= rt();
+  }
+};
+
+
+/* 
+   Gcalc_point_iterator simplifies the calculations on
+   the current slice of the Gcalc_scan_iterator.
+   One can walk through the points on the current slice.
+*/
+
+class Gcalc_point_iterator
+{
+protected:
+  const Gcalc_scan_iterator::point *sp;
+public:
+  Gcalc_point_iterator(const Gcalc_scan_iterator *scan_i):
+    sp(scan_i->get_b_slice())
+    {}
+
+  inline bool more() const { return sp != NULL; }
+  inline void operator++() { sp= sp->c_get_next(); }
+  inline const Gcalc_scan_iterator::point *point() const { return sp; }
+  inline const Gcalc_heap::Info *get_pi() const { return sp->pi; }
+  inline gcalc_shape_info get_shape() const { return sp->get_shape(); }
+  inline double get_x() const { return sp->x; }
+};
+
+#endif /*GCALC_SLICESCAN_INCLUDED*/
+
diff --git a/sql/gcalc_tools.cc b/sql/gcalc_tools.cc
new file mode 100644
index 00000000000..0e2970116ce
--- /dev/null
+++ b/sql/gcalc_tools.cc
@@ -0,0 +1,1156 @@
+/* Copyright (c) 2000, 2010 Oracle and/or its affiliates. All rights reserved.
+
+   This program is free software; you can redistribute it and/or modify
+   it under the terms of the GNU General Public License as published by
+   the Free Software Foundation; version 2 of the License.
+
+   This program is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+   GNU General Public License for more details.
+
+   You should have received a copy of the GNU General Public License
+   along with this program; if not, write to the Free Software
+   Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA */
+
+
+#include "mysql_priv.h"
+
+#ifdef HAVE_SPATIAL
+
+#include "gcalc_tools.h"
+#include "spatial.h"
+
+#define float_to_coord(d) ((double) d)
+
+
+/*
+  Adds new shape to the relation.
+  After that it can be used as an argument of an operation.
+*/
+
+gcalc_shape_info Gcalc_function::add_new_shape(uint32 shape_id,
+                                               shape_type shape_kind)
+{
+  shapes_buffer.q_append((uint32) shape_kind);
+  return n_shapes++;
+}
+
+
+/*
+  Adds new operation to the constructed relation.
+  To construct the complex relation one has to specify operations
+  in prefix style.
+*/
+
+void Gcalc_function::add_operation(op_type operation, uint32 n_operands)
+{
+  uint32 op_code= (uint32 ) operation + n_operands;
+  function_buffer.q_append(op_code);
+}
+
+
+/*
+  Sometimes the number of arguments is unknown at the moment the operation
+  is added. That allows to specify it later.
+*/
+
+void Gcalc_function::add_operands_to_op(uint32 operation_pos, uint32 n_operands)
+{
+  uint32 op_code= uint4korr(function_buffer.ptr() + operation_pos) + n_operands;
+  function_buffer.write_at_position(operation_pos, op_code);
+}
+
+
+/*
+  Just like the add_operation() but the result will be the inverted
+  value of an operation.
+*/
+
+void Gcalc_function::add_not_operation(op_type operation, uint32 n_operands)
+{
+  uint32 op_code= ((uint32) op_not | (uint32 ) operation) + n_operands;
+  function_buffer.q_append(op_code);
+}
+
+
+int Gcalc_function::single_shape_op(shape_type shape_kind, gcalc_shape_info *si)
+{
+  if (reserve_shape_buffer(1) || reserve_op_buffer(1))
+    return 1;
+  *si= add_new_shape(0, shape_kind);
+  add_operation(op_shape, *si);
+  return 0;
+}
+
+
+/*
+  Specify how many arguments we're going to have.
+*/
+
+int Gcalc_function::reserve_shape_buffer(uint n_shapes)
+{
+  return shapes_buffer.reserve(n_shapes * 4, 512);
+}
+
+
+/*
+  Specify how many operations we're going to have.
+*/
+
+int Gcalc_function::reserve_op_buffer(uint n_ops)
+{
+  return function_buffer.reserve(n_ops * 4, 512);
+}
+
+
+int Gcalc_function::alloc_states()
+{
+  if (function_buffer.reserve((n_shapes+1) * sizeof(int)))
+    return 1;
+  i_states= (int *) (function_buffer.ptr() + ALIGN_SIZE(function_buffer.length()));
+  return 0;
+}
+
+
+int Gcalc_function::count_internal()
+{
+  int c_op= uint4korr(cur_func);
+  op_type next_func= (op_type) (c_op & op_any);
+  int mask= (c_op & op_not) ? 1:0;
+  int n_ops= c_op & ~op_any;
+  int result;
+
+  cur_func+= 4;
+  if (next_func == op_shape)
+    return i_states[c_op & ~(op_any | op_not)] ^ mask;
+
+  result= count_internal();
+
+  while (--n_ops)
+  {
+    int next_res= count_internal();
+    switch (next_func)
+    {
+      case op_union:
+        result= result | next_res;
+        break;
+      case op_intersection:
+        result= result & next_res;
+        break;
+      case op_symdifference:
+        result= result ^ next_res;
+        break;
+      case op_difference:
+        result= result & !next_res;
+        break;
+      case op_backdifference:
+        result= !result & next_res;
+        break;
+      default:
+        DBUG_ASSERT(FALSE);
+    };
+  }
+
+  return result ^ mask;
+}
+
+
+/*
+  Clear the state of the object.
+*/
+
+void Gcalc_function::reset()
+{
+  n_shapes= 0;
+  shapes_buffer.length(0);
+  function_buffer.length(0);
+}
+
+
+int Gcalc_function::find_function(Gcalc_scan_iterator &scan_it)
+{
+  while (scan_it.more_points())
+  {
+    if (scan_it.step())
+      return -1;
+    Gcalc_scan_events ev= scan_it.get_event();
+    const Gcalc_scan_iterator::point *evpos= scan_it.get_event_position();
+    if (ev & (scev_point | scev_end | scev_two_ends))
+      continue;
+
+    clear_state();
+    for (Gcalc_point_iterator pit(&scan_it); pit.point() != evpos; ++pit)
+    {
+      gcalc_shape_info si= pit.point()->get_shape();
+      if ((get_shape_kind(si) == Gcalc_function::shape_polygon))
+        invert_state(si);
+    }
+    invert_state(evpos->get_shape());
+
+    if (ev == scev_intersection)
+    {
+      const Gcalc_scan_iterator::point *evnext= evpos->c_get_next();
+      if ((get_shape_kind(evpos->get_shape()) !=
+             Gcalc_function::shape_polygon)             ||
+          (get_shape_kind(evnext->get_shape()) !=
+             Gcalc_function::shape_polygon))
+        invert_state(evnext->get_shape());
+    }
+
+    if (count())
+      return 1;
+  }
+  return 0;
+}
+
+
+int Gcalc_operation_transporter::single_point(double x, double y)
+{
+  gcalc_shape_info si;
+  return m_fn->single_shape_op(Gcalc_function::shape_point, &si) ||
+         int_single_point(si, x, y);
+}
+
+
+int Gcalc_operation_transporter::start_line()
+{
+  int_start_line();
+  return m_fn->single_shape_op(Gcalc_function::shape_line, &m_si);
+}
+
+
+int Gcalc_operation_transporter::complete_line()
+{
+  int_complete_line();
+  return 0;
+}
+
+
+int Gcalc_operation_transporter::start_poly()
+{
+  int_start_poly();
+  return m_fn->single_shape_op(Gcalc_function::shape_polygon, &m_si);
+}
+
+
+int Gcalc_operation_transporter::complete_poly()
+{
+  int_complete_poly();
+  return 0;
+}
+
+
+int Gcalc_operation_transporter::start_ring()
+{
+  int_start_ring();
+  return 0;
+}
+
+
+int Gcalc_operation_transporter::complete_ring()
+{
+  int_complete_ring();
+  return 0;
+}
+
+
+int Gcalc_operation_transporter::add_point(double x, double y)
+{
+  return int_add_point(m_si, x, y);
+}
+
+
+int Gcalc_operation_transporter::start_collection(int n_objects)
+{
+  if (m_fn->reserve_shape_buffer(n_objects) || m_fn->reserve_op_buffer(1))
+        return 1;
+  m_fn->add_operation(Gcalc_function::op_union, n_objects);
+  return 0;
+}
+
+
+int Gcalc_result_receiver::start_shape(Gcalc_function::shape_type shape)
+{
+  if (buffer.reserve(4*2, 512))
+    return 1;
+  cur_shape= shape;
+  shape_pos= buffer.length();
+  buffer.length(shape_pos + ((shape == Gcalc_function::shape_point) ? 4:8));
+  n_points= 0;
+  shape_area= 0.0;
+
+  return 0;
+}
+
+
+int Gcalc_result_receiver::add_point(double x, double y)
+{
+  if (n_points && x == prev_x && y == prev_y)
+    return 0;
+
+  if (!n_points++)
+  {
+    prev_x= first_x= x;
+    prev_y= first_y= y;
+    return 0;
+  }
+
+  shape_area+= prev_x*y - prev_y*x;
+
+  if (buffer.reserve(8*2, 512))
+    return 1;
+  buffer.q_append(prev_x);
+  buffer.q_append(prev_y);
+  prev_x= x;
+  prev_y= y;
+  return 0;
+}
+
+
+int Gcalc_result_receiver::complete_shape()
+{
+  if (n_points == 0)
+  {
+    buffer.length(shape_pos);
+    return 0;
+  }
+  if (n_points == 1)
+  {
+    if (cur_shape != Gcalc_function::shape_point)
+    {
+      cur_shape= Gcalc_function::shape_point;
+      buffer.length(buffer.length()-4);
+    }
+  }
+  else
+  {
+    DBUG_ASSERT(cur_shape != Gcalc_function::shape_point);
+    if (cur_shape == Gcalc_function::shape_hole)
+    {
+      shape_area+= prev_x*first_y - prev_y*first_x;
+      if (fabs(shape_area) < 1e-8)
+      {
+        buffer.length(shape_pos);
+        return 0;
+      }
+    }
+
+    if ((cur_shape == Gcalc_function::shape_polygon ||
+          cur_shape == Gcalc_function::shape_hole) &&
+        prev_x == first_x && prev_y == first_y)
+    {
+      n_points--;
+      buffer.write_at_position(shape_pos+4, n_points);
+      goto do_complete;
+    }
+    buffer.write_at_position(shape_pos+4, n_points);
+  }
+
+  if (buffer.reserve(8*2, 512))
+    return 1;
+  buffer.q_append(prev_x);
+  buffer.q_append(prev_y);
+  
+do_complete:
+  buffer.write_at_position(shape_pos, (uint32) cur_shape);
+
+  if (!n_shapes++)
+  {
+    DBUG_ASSERT(cur_shape != Gcalc_function::shape_hole);
+    common_shapetype= cur_shape;
+  }
+  else if (cur_shape == Gcalc_function::shape_hole)
+  {
+    ++n_holes;
+  }
+  else if (!collection_result && (cur_shape != common_shapetype))
+  {
+      collection_result= true;
+  }
+  return 0;
+}
+
+
+int Gcalc_result_receiver::single_point(double x, double y)
+{
+  return start_shape(Gcalc_function::shape_point) ||
+         add_point(x, y) ||
+         complete_shape();
+}
+
+
+int Gcalc_result_receiver::done()
+{
+  return 0;
+}
+
+
+void Gcalc_result_receiver::reset()
+{
+  buffer.length(0);
+  collection_result= FALSE;
+  n_shapes= n_holes= 0;
+}
+
+
+int Gcalc_result_receiver::get_result_typeid()
+{
+  if (!n_shapes)
+    return 0;
+
+  if (collection_result)
+    return Geometry::wkb_geometrycollection;
+  switch (common_shapetype)
+  {
+    case Gcalc_function::shape_polygon:
+      return (n_shapes - n_holes == 1) ?
+              Geometry::wkb_polygon : Geometry::wkb_multipolygon;
+    case Gcalc_function::shape_point:
+      return (n_shapes == 1) ? Geometry::wkb_point : Geometry::wkb_multipoint;
+    case Gcalc_function::shape_line:
+      return (n_shapes == 1) ? Geometry::wkb_linestring :
+                               Geometry::wkb_multilinestring;
+    default:
+      DBUG_ASSERT(0);
+  }
+  return 0;
+}
+
+
+int Gcalc_result_receiver::move_hole(uint32 dest_position, uint32 source_position,
+                                     uint32 *new_dest_position)
+{
+  char *ptr;
+  int source_len;
+  if (dest_position == source_position)
+  {
+    *new_dest_position= position();
+    return 0;
+  }
+
+  source_len= buffer.length() - source_position;
+  if (buffer.reserve(source_len, MY_ALIGN(source_len, 512)))
+    return 1;
+
+  ptr= (char *) buffer.ptr();
+  memmove(ptr + dest_position + source_len, ptr + dest_position,
+          buffer.length() - dest_position);
+  memcpy(ptr + dest_position, ptr + buffer.length(), source_len);
+  *new_dest_position= dest_position + source_len;
+  return 0;
+}
+
+
+Gcalc_operation_reducer::Gcalc_operation_reducer(size_t blk_size) :
+  Gcalc_dyn_list(blk_size, sizeof(res_point)),
+  m_res_hook((Gcalc_dyn_list::Item **)&m_result),
+  m_first_active_thread(NULL)
+{}
+
+
+void Gcalc_operation_reducer::init(Gcalc_function *fn, modes mode)
+{
+  m_fn= fn;
+  m_mode= mode;
+  m_first_active_thread= NULL;
+}
+
+
+Gcalc_operation_reducer::
+Gcalc_operation_reducer(Gcalc_function *fn, modes mode, size_t blk_size) :
+  Gcalc_dyn_list(blk_size, sizeof(res_point)),
+  m_res_hook((Gcalc_dyn_list::Item **)&m_result)
+{
+  init(fn, mode);
+}
+
+
+inline int Gcalc_operation_reducer::continue_range(active_thread *t,
+						const Gcalc_heap::Info *p)
+{
+  DBUG_ASSERT(t->result_range);
+  res_point *rp= add_res_point();
+  if (!rp)
+    return 1;
+  rp->glue= NULL;
+  rp->down= t->rp;
+  t->rp->up= rp;
+  rp->intersection_point= false;
+  rp->pi= p;
+  t->rp= rp;
+  return 0;
+}
+
+
+inline int Gcalc_operation_reducer::continue_i_range(active_thread *t,
+						  const Gcalc_heap::Info *p,
+						  double x, double y)
+{
+  DBUG_ASSERT(t->result_range);
+  res_point *rp= add_res_point();
+  if (!rp)
+    return 1;
+  rp->glue= NULL;
+  rp->down= t->rp;
+  t->rp->up= rp;
+  rp->intersection_point= true;
+  rp->x= x;
+  rp->pi= p;
+  rp->y= y;
+  t->rp= rp;
+  return 0;
+}
+
+inline int Gcalc_operation_reducer::start_range(active_thread *t,
+					     const Gcalc_heap::Info *p)
+{
+  res_point *rp= add_res_point();
+  if (!rp)
+    return 1;
+  rp->glue= rp->down= NULL;
+  rp->intersection_point= false;
+  rp->pi= p;
+  t->result_range= 1;
+  t->rp= rp;
+  return 0;
+}
+
+inline int Gcalc_operation_reducer::start_i_range(active_thread *t,
+					       const Gcalc_heap::Info *p,
+					       double x, double y)
+{
+  res_point *rp= add_res_point();
+  if (!rp)
+    return 1;
+  rp->glue= rp->down= NULL;
+  rp->intersection_point= true;
+  rp->x= x;
+  rp->y= y;
+  rp->pi= p;
+  t->result_range= 1;
+  t->rp= rp;
+  return 0;
+}
+
+inline int Gcalc_operation_reducer::end_range(active_thread *t,
+					   const Gcalc_heap::Info *p)
+{
+  res_point *rp= add_res_point();
+  if (!rp)
+    return 1;
+  rp->glue= rp->up= NULL;
+  rp->down= t->rp;
+  rp->intersection_point= false;
+  rp->pi= p;
+  t->rp->up= rp;
+  t->result_range= 0;
+  return 0;
+}
+
+inline int Gcalc_operation_reducer::end_i_range(active_thread *t,
+					     const Gcalc_heap::Info *p,
+					     double x, double y)
+{
+  res_point *rp= add_res_point();
+  if (!rp)
+    return 1;
+  rp->glue= rp->up= NULL;
+  rp->down= t->rp;
+  rp->intersection_point= true;
+  rp->x= x;
+  rp->pi= p;
+  rp->y= y;
+  t->rp->up= rp;
+  t->result_range= 0;
+  return 0;
+}
+
+int Gcalc_operation_reducer::start_couple(active_thread *t0, active_thread *t1,
+				       const Gcalc_heap::Info *p,
+                                       const active_thread *prev_range)
+{
+  res_point *rp0, *rp1;
+  if (!(rp0= add_res_point()) || !(rp1= add_res_point()))
+    return 1;
+  rp0->glue= rp1;
+  rp1->glue= rp0;
+  rp0->intersection_point= rp1->intersection_point= false;
+  rp0->down= rp1->down= NULL;
+  rp0->pi= rp1->pi= p;
+  t0->rp= rp0;
+  t1->rp= rp1;
+  if (prev_range)
+  {
+    rp0->outer_poly= prev_range->thread_start;
+    t1->thread_start= prev_range->thread_start;
+  }
+  else
+  {
+    rp0->outer_poly= 0;
+    t0->thread_start= rp0;
+  }
+  return 0;
+}
+
+int Gcalc_operation_reducer::start_i_couple(active_thread *t0, active_thread *t1,
+					 const Gcalc_heap::Info *p0,
+					 const Gcalc_heap::Info *p1,
+					 double x, double y,
+                                         const active_thread *prev_range)
+{
+  res_point *rp0, *rp1;
+  if (!(rp0= add_res_point()) || !(rp1= add_res_point()))
+    return 1;
+  rp0->glue= rp1;
+  rp1->glue= rp0;
+  rp0->pi= p0;
+  rp1->pi= p1;
+  rp0->intersection_point= rp1->intersection_point= true;
+  rp0->down= rp1->down= NULL;
+  rp0->x= rp1->x= x;
+  rp0->y= rp1->y= y;
+  t0->result_range= t1->result_range= 1;
+  t0->rp= rp0;
+  t1->rp= rp1;
+  if (prev_range)
+  {
+    rp0->outer_poly= prev_range->thread_start;
+    t1->thread_start= prev_range->thread_start;
+  }
+  else
+  {
+    rp0->outer_poly= 0;
+    t0->thread_start= rp0;
+  }
+  return 0;
+}
+
+int Gcalc_operation_reducer::end_couple(active_thread *t0, active_thread *t1,
+				     const Gcalc_heap::Info *p)
+{
+  res_point *rp0, *rp1;
+  DBUG_ASSERT(t1->result_range);
+  if (!(rp0= add_res_point()) || !(rp1= add_res_point()))
+    return 1;
+  rp0->down= t0->rp;
+  rp1->down= t1->rp;
+  rp1->glue= rp0;
+  rp0->glue= rp1;
+  rp0->up= rp1->up= NULL;
+  t0->rp->up= rp0;
+  t1->rp->up= rp1;
+  rp0->intersection_point= rp1->intersection_point= false;
+  rp0->pi= rp1->pi= p;
+  t0->result_range= t1->result_range= 0;
+  return 0;
+}
+
+int Gcalc_operation_reducer::end_i_couple(active_thread *t0, active_thread *t1,
+				       const Gcalc_heap::Info *p0,
+				       const Gcalc_heap::Info *p1,
+				       double x, double y)
+{
+  res_point *rp0, *rp1;
+  if (!(rp0= add_res_point()) || !(rp1= add_res_point()))
+    return 1;
+  rp0->down= t0->rp;
+  rp1->down= t1->rp;
+  rp0->pi= p0;
+  rp1->pi= p1;
+  rp1->glue= rp0;
+  rp0->glue= rp1;
+  rp0->up= rp1->up= NULL;
+  rp0->intersection_point= rp1->intersection_point= true;
+  rp0->x= rp1->x= x;
+  rp0->y= rp1->y= y;
+  t0->result_range= t1->result_range= 0;
+  t0->rp->up= rp0;
+  t1->rp->up= rp1;
+  return 0;
+}
+
+int Gcalc_operation_reducer::add_single_point(const Gcalc_heap::Info *p)
+{
+  res_point *rp= add_res_point();
+  if (!rp)
+    return 1;
+  rp->glue= rp->up= rp->down= NULL;
+  rp->intersection_point= false;
+  rp->pi= p;
+  rp->x= p->x;
+  rp->y= p->y;
+  return 0;
+}
+
+int Gcalc_operation_reducer::add_i_single_point(const Gcalc_heap::Info *p,
+					     double x, double y)
+{
+  res_point *rp= add_res_point();
+  if (!rp)
+    return 1;
+  rp->glue= rp->up= rp->down= NULL;
+  rp->intersection_point= true;
+  rp->x= x;
+  rp->pi= p;
+  rp->y= y;
+  return 0;
+}
+
+int Gcalc_operation_reducer::
+handle_lines_intersection(active_thread *t0, active_thread *t1,
+			  const Gcalc_heap::Info *p0, const Gcalc_heap::Info *p1,
+			  double x, double y)
+{
+  m_fn->invert_state(p0->shape);
+  m_fn->invert_state(p1->shape);
+  int intersection_state= m_fn->count();
+  if ((t0->result_range | t1->result_range) == intersection_state)
+    return 0;
+
+  if (t0->result_range &&
+      (end_i_range(t0, p1, x, y) || start_i_range(t0, p1, x, y)))
+    return 1;
+
+  if (t1->result_range &&
+      (end_i_range(t1, p0, x, y) || start_i_range(t1, p0, x, y)))
+    return 1;
+
+  if (intersection_state &&
+      add_i_single_point(p0, x, y))
+    return 1;
+    
+  return 0;
+}
+
+inline int Gcalc_operation_reducer::
+handle_line_polygon_intersection(active_thread *l, const Gcalc_heap::Info *pl,
+				 int line_state, int poly_state,
+				 double x, double y)
+{
+  int range_after= ~poly_state & line_state;
+  if (l->result_range == range_after)
+    return 0;
+  return range_after ? start_i_range(l, pl, x, y) : end_i_range(l, pl, x, y);
+}
+
+static inline void switch_athreads(Gcalc_operation_reducer::active_thread *t0,
+				   Gcalc_operation_reducer::active_thread *t1,
+				   Gcalc_dyn_list::Item **hook)
+{
+  *hook= t1;
+  t0->next= t1->next;
+  t1->next= t0;
+}
+
+inline int Gcalc_operation_reducer::
+handle_polygons_intersection(active_thread *t0, active_thread *t1,
+			     Gcalc_dyn_list::Item **t_hook,
+			     const Gcalc_heap::Info *p0,
+			     const Gcalc_heap::Info *p1,
+			     int prev_state, double x, double y,
+                             const active_thread *prev_range)
+{
+  m_fn->invert_state(p0->shape);
+  int state_11= m_fn->count();
+  m_fn->invert_state(p1->shape);
+  int state_2= m_fn->count();
+  int state_01= prev_state ^ t0->result_range;
+  if ((prev_state == state_01) && (prev_state == state_2))
+  {
+    if (state_11 == prev_state)
+    {
+      switch_athreads(t0, t1, t_hook);
+      return 0;
+    }
+    return start_i_couple(t0, t1, p0, p1, x, y, prev_range);
+  }
+  if (prev_state == state_2)
+  {
+    if (state_01 == state_11)
+    {
+      if (m_mode & polygon_selfintersections_allowed)
+      {
+	switch_athreads(t0, t1, t_hook);
+	return 0;
+      }
+      if (prev_state != (m_mode & prefer_big_with_holes))
+	return continue_i_range(t0, p0, x, y) || continue_i_range(t1, p1, x, y);
+      return end_i_couple(t0, t1, p0, p1, x, y) ||
+	start_i_couple(t0, t1, p0, p1, x, y, prev_range);
+    }
+    else
+      return end_i_couple(t0, t1, p0, p1, x, y);
+  }
+  if (state_01 ^ state_11)
+  {
+    switch_athreads(t0, t1, t_hook);
+    return 0;
+  }
+
+  active_thread *thread_to_continue;
+  const Gcalc_heap::Info *way_to_go;
+  if (prev_state == state_01)
+  {
+    thread_to_continue= t1;
+    way_to_go= p1;
+  }
+  else 
+  {
+    thread_to_continue= t0;
+    way_to_go= p0;
+  }
+  return continue_i_range(thread_to_continue, way_to_go, x, y);
+}
+
+int Gcalc_operation_reducer::count_slice(Gcalc_scan_iterator *si)
+{
+  Gcalc_point_iterator pi(si);
+  active_thread *cur_t= m_first_active_thread;
+  Gcalc_dyn_list::Item **at_hook= (Gcalc_dyn_list::Item **)&m_first_active_thread;
+  const active_thread *prev_range;
+  int prev_state;
+
+  if (si->get_event() & (scev_point | scev_end | scev_two_ends))
+  {
+    for (; pi.point() != si->get_event_position(); ++pi, cur_t= cur_t->get_next())
+      at_hook= &cur_t->next;
+
+    switch (si->get_event())
+    {
+      case scev_point:
+      {
+        if (cur_t->result_range &&
+            continue_range(cur_t, pi.get_pi()))
+          return 1;
+        break;
+      }
+      case scev_end:
+      {
+        if (cur_t->result_range &&
+            end_range(cur_t, pi.get_pi()))
+          return 1;
+        *at_hook= cur_t->next;
+        free_item(cur_t);
+        break;
+      }
+      case scev_two_ends:
+      {
+        active_thread *cur_t1= cur_t->get_next();
+        if (cur_t->result_range &&
+            end_couple(cur_t, cur_t1, pi.get_pi()))
+          return 1;
+
+        *at_hook= cur_t1->next;
+        free_list(cur_t, &cur_t1->next);
+        break;
+      }
+      default:
+        DBUG_ASSERT(0);
+    }
+    return 0;
+  }
+
+  prev_state= 0;
+  prev_range= 0;
+
+  m_fn->clear_state();
+  for (; pi.point() != si->get_event_position(); ++pi, cur_t= cur_t->get_next())
+  {
+    if (m_fn->get_shape_kind(pi.get_shape()) == Gcalc_function::shape_polygon)
+    {
+      m_fn->invert_state(pi.get_shape());
+      prev_state^= cur_t->result_range;
+    }
+    at_hook= &cur_t->next;
+    if (cur_t->result_range)
+      prev_range= prev_state ? cur_t : 0;
+  }
+
+  switch (si->get_event())
+  {
+  case scev_thread:
+  {
+    active_thread *new_t= new_active_thread();
+    if (!new_t)
+      return 1;
+    m_fn->invert_state(pi.get_shape());
+    new_t->result_range= prev_state ^ m_fn->count();
+    new_t->next= *at_hook;
+    *at_hook= new_t;
+    if (new_t->result_range &&
+	start_range(new_t, pi.get_pi()))
+      return 1;
+    break;
+  }
+  case scev_two_threads:
+  {
+    active_thread *new_t0, *new_t1;
+    int fn_result;
+    if (!(new_t0= new_active_thread()) || !(new_t1= new_active_thread()))
+      return 1;
+    
+    m_fn->invert_state(pi.get_shape());
+    fn_result= m_fn->count();
+    new_t0->result_range= new_t1->result_range= prev_state ^ fn_result;
+    new_t1->next= *at_hook;
+    new_t0->next= new_t1;
+    *at_hook= new_t0;
+    if (new_t0->result_range &&
+	start_couple(new_t0, new_t1, pi.get_pi(), prev_range))
+      return 1;
+    break;
+  }
+  case scev_intersection:
+  {
+    active_thread *cur_t1= cur_t->get_next();
+    const Gcalc_heap::Info *p0, *p1;
+    p0= pi.get_pi();
+    ++pi;
+    p1= pi.get_pi();
+    bool line0= m_fn->get_shape_kind(p0->shape) == Gcalc_function::shape_line;
+    bool line1= m_fn->get_shape_kind(p1->shape) == Gcalc_function::shape_line;
+
+    if (!line0 && !line1) /* two polygons*/
+    {
+      if (handle_polygons_intersection(cur_t, cur_t1, at_hook, p0, p1,
+				       prev_state, pi.get_x(), si->get_y(),
+                                       prev_range))
+	return 1;
+    }
+    else if (line0 && line1)
+    {
+      if (!prev_state &&
+	  handle_lines_intersection(cur_t, cur_t1,
+				    p0, p1, pi.get_x(), si->get_y()))
+	return 1;
+      switch_athreads(cur_t, cur_t1, at_hook);
+    }
+    else
+    {
+      int poly_state;
+      int line_state;
+      const Gcalc_heap::Info *line;
+      active_thread *line_t;
+      m_fn->invert_state(p0->shape);
+      if (line0)
+      {
+	line_state= m_fn->count();
+	poly_state= prev_state;
+	line= p0;
+	line_t= cur_t1;
+      }
+      else
+      {
+	poly_state= m_fn->count();
+	m_fn->invert_state(p1->shape);
+	line_state= m_fn->count();
+	line= p1;
+	line_t= cur_t;
+      }
+      if (handle_line_polygon_intersection(line_t, line,
+					   line_state, poly_state,
+					   pi.get_x(), si->get_y()))
+	return 1;
+      switch_athreads(cur_t, cur_t1, at_hook);
+    }
+    break;
+  }
+  case scev_single_point:
+  {
+    m_fn->invert_state(pi.get_shape());
+    if ((prev_state ^ m_fn->count()) &&
+	add_single_point(pi.get_pi()))
+      return 1;
+    break;
+  }
+  default:
+    DBUG_ASSERT(0);
+  }
+
+  return 0;
+}
+
+int Gcalc_operation_reducer::count_all(Gcalc_heap *hp)
+{
+  Gcalc_scan_iterator si;
+  si.init(hp);
+  while (si.more_points())
+  {
+    if (si.step())
+      return 1;
+    if (count_slice(&si))
+      return 1;
+  }
+  return 0;
+}
+
+inline void Gcalc_operation_reducer::free_result(res_point *res)
+{
+  if ((*res->prev_hook= res->next))
+  {
+    res->get_next()->prev_hook= res->prev_hook;
+  }
+  free_item(res);
+}
+
+
+inline int Gcalc_operation_reducer::get_single_result(res_point *res,
+						   Gcalc_result_receiver *storage)
+{
+  if (res->intersection_point)
+  {
+    if (storage->single_point(float_to_coord(res->x),
+			      float_to_coord(res->y)))
+      return 1;
+  }
+  else
+    if (storage->single_point(res->x, res->y))
+      return 1;
+  free_result(res);
+  return 0;
+}
+
+
+int Gcalc_operation_reducer::get_result_thread(res_point *cur,
+                                               Gcalc_result_receiver *storage,
+                                               int move_upward)
+{
+  res_point *next;
+  bool glue_step= false;
+  res_point *first_poly_node= cur;
+  double x, y;
+  while (cur)
+  {
+    if (!glue_step)
+    {
+      if (cur->intersection_point)
+      {
+        x= float_to_coord(cur->x);
+        y= float_to_coord(cur->y);
+      }
+      else
+      {
+	x= cur->pi->x;
+        y= cur->pi->y;
+      }
+      if (storage->add_point(x, y))
+        return 1;
+    }
+    
+    next= move_upward ? cur->up : cur->down;
+    if (!next && !glue_step)
+    {
+      next= cur->glue;
+      move_upward^= 1;
+      glue_step= true;
+      if (next)
+	next->glue= NULL;
+    }
+    else
+      glue_step= false;
+
+    cur->first_poly_node= first_poly_node;
+    free_result(cur);
+    cur= next;
+  }
+  return 0;
+}
+
+
+int Gcalc_operation_reducer::get_polygon_result(res_point *cur,
+                                                Gcalc_result_receiver *storage)
+{
+  res_point *glue= cur->glue;
+  glue->up->down= NULL;
+  free_result(glue);
+  return get_result_thread(cur, storage, 1) ||
+         storage->complete_shape();
+}
+
+
+int Gcalc_operation_reducer::get_line_result(res_point *cur,
+                                             Gcalc_result_receiver *storage)
+{
+  res_point *next;
+  int move_upward= 1;
+  if (cur->glue)
+  {
+    /* Here we have to find the beginning of the line */
+    next= cur->up;
+    move_upward= 1;
+    while (next)
+    {
+      cur= next;
+      next= move_upward ? next->up : next->down;
+      if (!next)
+      {
+	next= cur->glue;
+	move_upward^= 1;
+      }
+    }
+  }
+
+  return get_result_thread(cur, storage, move_upward) ||
+         storage->complete_shape();
+}
+
+
+int Gcalc_operation_reducer::get_result(Gcalc_result_receiver *storage)
+{
+  *m_res_hook= NULL;
+  while (m_result)
+  {
+    if (!m_result->up)
+    {
+      if (get_single_result(m_result, storage))
+	return 1;
+      continue;
+    }
+    Gcalc_function::shape_type shape= m_fn->get_shape_kind(m_result->pi->shape);
+    if (shape == Gcalc_function::shape_polygon)
+    {
+      if (m_result->outer_poly)
+      {
+        uint32 *insert_position, hole_position;
+        insert_position= &m_result->outer_poly->first_poly_node->poly_position;
+        DBUG_ASSERT(*insert_position);
+        hole_position= storage->position();
+        storage->start_shape(Gcalc_function::shape_hole);
+        if (get_polygon_result(m_result, storage) ||
+            storage->move_hole(*insert_position, hole_position,
+                               insert_position))
+          return 1;
+      }
+      else
+      {
+        uint32 *poly_position= &m_result->poly_position;
+        storage->start_shape(Gcalc_function::shape_polygon);
+        if (get_polygon_result(m_result, storage))
+          return 1;
+        *poly_position= storage->position();
+      }
+    }
+    else
+    {
+      storage->start_shape(shape);
+      if (get_line_result(m_result, storage))
+	return 1;
+    }
+  }
+  
+  m_res_hook= (Gcalc_dyn_list::Item **)&m_result;
+  storage->done();
+  return 0;
+}
+
+
+void Gcalc_operation_reducer::reset()
+{
+  free_list(m_result, m_res_hook);
+  m_res_hook= (Gcalc_dyn_list::Item **)&m_result;
+  free_list(m_first_active_thread);
+}
+
+#endif /*HAVE_SPATIAL*/
+
diff --git a/sql/gcalc_tools.h b/sql/gcalc_tools.h
new file mode 100644
index 00000000000..7df236618fb
--- /dev/null
+++ b/sql/gcalc_tools.h
@@ -0,0 +1,306 @@
+/* Copyright (c) 2000, 2010 Oracle and/or its affiliates. All rights reserved.
+
+   This program is free software; you can redistribute it and/or modify
+   it under the terms of the GNU General Public License as published by
+   the Free Software Foundation; version 2 of the License.
+
+   This program is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+   GNU General Public License for more details.
+
+   You should have received a copy of the GNU General Public License
+   along with this program; if not, write to the Free Software
+   Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA */
+
+
+#ifndef GCALC_TOOLS_INCLUDED
+#define GCALC_TOOLS_INCLUDED
+
+#include "gcalc_slicescan.h"
+
+
+/*
+  The Gcalc_function class objects are used to check for a binary relation.
+  The relation can be constructed with the prefix notation using predicates as
+        op_not (as !A)
+        op_union ( A || B || C... )
+        op_intersection ( A && B && C ... )
+        op_symdifference ( A+B+C+... == 1 )
+        op_difference ( A && !(B||C||..))
+  with the calls of the add_operation(operation, n_operands) method.
+  The relation is calculated over a set of shapes, that in turn have
+  to be added with the add_new_shape() method. All the 'shapes' can
+  be set to 0 with clear_shapes() method and single value
+  can be changed with the invert_state() method.
+  Then the value of the relation can be calculated with the count() method.
+  Frequently used method is find_function(Gcalc_scan_iterator it) that
+  iterates through the 'it' until the relation becomes TRUE.
+*/
+
+class Gcalc_function
+{
+private:
+  String shapes_buffer;
+  String function_buffer;
+  const char *cur_func;
+  int *i_states;
+  uint32 cur_object_id;
+  uint n_shapes;
+  int count_internal();
+public:
+  enum op_type
+  {
+    op_shape= 0,
+    op_not= 0x80000000,
+    op_union= 0x10000000,
+    op_intersection= 0x20000000,
+    op_symdifference= 0x30000000,
+    op_difference= 0x40000000,
+    op_backdifference= 0x50000000,
+    op_any= 0x70000000
+  };
+  enum shape_type
+  {
+    shape_point= 0,
+    shape_line= 1,
+    shape_polygon= 2,
+    shape_hole= 3
+  };
+  Gcalc_function() : n_shapes(0) {}
+  gcalc_shape_info add_new_shape(uint32 shape_id, shape_type shape_kind);
+  /*
+    Adds the leaf operation that returns the shape value.
+    Also adds the shape to the list of operands.
+  */
+  int single_shape_op(shape_type shape_kind, gcalc_shape_info *si);
+  void add_operation(op_type operation, uint32 n_operands);
+  void add_not_operation(op_type operation, uint32 n_operands);
+  uint32 get_next_operation_pos() { return function_buffer.length(); }
+  void add_operands_to_op(uint32 operation_pos, uint32 n_operands);
+  void set_cur_obj(uint32 cur_obj) { cur_object_id= cur_obj; }
+  int reserve_shape_buffer(uint n_shapes);
+  int reserve_op_buffer(uint n_ops);
+  uint get_nshapes() const { return n_shapes; }
+  shape_type get_shape_kind(gcalc_shape_info si) const
+  {
+    return (shape_type) uint4korr(shapes_buffer.ptr() + (si*4));
+  }
+
+  void set_states(int *shape_states) { i_states= shape_states; }
+  int alloc_states();
+  void invert_state(gcalc_shape_info shape) { i_states[shape]^= 1; }
+  int get_state(gcalc_shape_info shape) { return i_states[shape]; }
+  int count()
+  {
+    cur_func= function_buffer.ptr();
+    return count_internal();
+  }
+  void clear_state() { bzero(i_states, n_shapes * sizeof(int)); }
+  void reset();
+
+  int find_function(Gcalc_scan_iterator &scan_it);
+};
+
+
+/*
+  Gcalc_operation_transporter class extends the Gcalc_shape_transporter.
+  In addition to the parent's functionality, it fills the Gcalc_function
+  object so it has the function that determines the proper shape.
+  For example Multipolyline will be represented as an union of polylines.
+*/
+
+class Gcalc_operation_transporter : public Gcalc_shape_transporter
+{
+protected:
+  Gcalc_function *m_fn;
+  gcalc_shape_info m_si;
+public:
+  Gcalc_operation_transporter(Gcalc_function *fn, Gcalc_heap *heap) :
+    Gcalc_shape_transporter(heap), m_fn(fn) {}
+
+  int single_point(double x, double y);
+  int start_line();
+  int complete_line();
+  int start_poly();
+  int complete_poly();
+  int start_ring();
+  int complete_ring();
+  int add_point(double x, double y);
+  int start_collection(int n_objects);
+};
+
+
+/*
+   When we calculate the result of an spatial operation like
+   Union or Intersection, we receive vertexes of the result
+   one-by-one, and probably need to treat them in variative ways.
+   So, the Gcalc_result_receiver class designed to get these
+   vertexes and construct shapes/objects out of them.
+   and to store the result in an appropriate format
+*/
+
+class Gcalc_result_receiver
+{
+  String buffer;
+  uint32 n_points;
+  Gcalc_function::shape_type common_shapetype;
+  bool collection_result;
+  uint32 n_shapes;
+  uint32 n_holes;
+
+  Gcalc_function::shape_type cur_shape;
+  uint32 shape_pos;
+  double first_x, first_y, prev_x, prev_y;
+  double shape_area;
+public:
+  Gcalc_result_receiver() : collection_result(FALSE), n_shapes(0), n_holes(0)
+    {}
+  int start_shape(Gcalc_function::shape_type shape);
+  int add_point(double x, double y);
+  int complete_shape();
+  int single_point(double x, double y);
+  int done();
+  void reset();
+
+  const char *result() { return buffer.ptr(); }
+  uint length() { return buffer.length(); }
+  int get_nshapes() { return n_shapes; }
+  int get_nholes() { return n_holes; }
+  int get_result_typeid();
+  uint32 position() { return buffer.length(); }
+  int move_hole(uint32 dest_position, uint32 source_position,
+                uint32 *new_dest_position);
+};
+
+
+/*
+  Gcalc_operation_reducer class incapsulates the spatial
+  operation functionality. It analyses the slices generated by
+  the slicescan and calculates the shape of the result defined
+  by some Gcalc_function.
+*/
+
+class Gcalc_operation_reducer : public Gcalc_dyn_list
+{
+public:
+  enum modes
+  {
+    /* Numeric values important here - careful with changing */
+    default_mode= 0,
+    prefer_big_with_holes= 1,
+    polygon_selfintersections_allowed= 2,  /* allowed in the result */
+    line_selfintersections_allowed= 4      /* allowed in the result */
+  };
+
+  Gcalc_operation_reducer(size_t blk_size=8192);
+  void init(Gcalc_function *fn, modes mode= default_mode);
+  Gcalc_operation_reducer(Gcalc_function *fn, modes mode= default_mode,
+		       size_t blk_size=8192);
+  int count_slice(Gcalc_scan_iterator *si);
+  int count_all(Gcalc_heap *hp);
+  int get_result(Gcalc_result_receiver *storage);
+  void reset();
+
+  class res_point : public Gcalc_dyn_list::Item
+  {
+  public:
+    bool intersection_point;
+    double x,y;
+    res_point *up;
+    res_point *down;
+    res_point *glue;
+    union
+    {
+      const Gcalc_heap::Info *pi;
+      res_point *first_poly_node;
+    };
+    union
+    {
+      res_point *outer_poly;
+      uint32 poly_position;
+    };
+    Gcalc_dyn_list::Item **prev_hook;
+    res_point *get_next() { return (res_point *)next; }
+  };
+
+  class active_thread : public Gcalc_dyn_list::Item
+  {
+  public:
+    res_point *rp;
+    int result_range;
+    res_point *thread_start;
+    active_thread *get_next() { return (active_thread *)next; }
+  };
+
+protected:
+  Gcalc_function *m_fn;
+  Gcalc_dyn_list::Item **m_res_hook;
+  res_point *m_result;
+  int m_mode;
+
+  res_point *result_heap;
+  active_thread *m_first_active_thread;
+
+  res_point *add_res_point()
+  {
+    res_point *result= (res_point *)new_item();
+    *m_res_hook= result;
+    result->prev_hook= m_res_hook;
+    m_res_hook= &result->next;
+    return result;
+  }
+
+  active_thread *new_active_thread() { return (active_thread *)new_item(); }
+
+private:
+  int continue_range(active_thread *t, const Gcalc_heap::Info *p);
+  int continue_i_range(active_thread *t, const Gcalc_heap::Info *p,
+		       double x, double y);
+  int start_range(active_thread *t, const Gcalc_heap::Info *p);
+  int start_i_range(active_thread *t, const Gcalc_heap::Info *p,
+		    double x, double y);
+  int end_range(active_thread *t, const Gcalc_heap::Info *p);
+  int end_i_range(active_thread *t, const Gcalc_heap::Info *p,
+		  double x, double y);
+  int start_couple(active_thread *t0, active_thread *t1,const Gcalc_heap::Info *p,
+                     const active_thread *prev_range);
+  int start_i_couple(active_thread *t0, active_thread *t1,
+		     const Gcalc_heap::Info *p0,
+		     const Gcalc_heap::Info *p1,
+		     double x, double y,
+                     const active_thread *prev_range);
+  int end_couple(active_thread *t0, active_thread *t1, const Gcalc_heap::Info *p);
+  int end_i_couple(active_thread *t0, active_thread *t1,
+		   const Gcalc_heap::Info *p0,
+		   const Gcalc_heap::Info *p1,
+		   double x, double y);
+  int add_single_point(const Gcalc_heap::Info *p);
+  int add_i_single_point(const Gcalc_heap::Info *p, double x, double y);
+
+  int handle_lines_intersection(active_thread *t0, active_thread *t1,
+				const Gcalc_heap::Info *p0,
+				const Gcalc_heap::Info *p1,
+				double x, double y);
+  int handle_polygons_intersection(active_thread *t0, active_thread *t1,
+				   Gcalc_dyn_list::Item **t_hook,
+				   const Gcalc_heap::Info *p0,
+				   const Gcalc_heap::Info *p1,
+				   int prev_state, double x, double y,
+                                   const active_thread *prev_range);
+  int handle_line_polygon_intersection(active_thread *l,
+				       const Gcalc_heap::Info *pl,
+				       int line_state, int poly_state,
+				       double x, double y);
+
+  int get_single_result(res_point *res, Gcalc_result_receiver *storage);
+  int get_result_thread(res_point *cur, Gcalc_result_receiver *storage,
+			int move_upward);
+  int get_polygon_result(res_point *cur, Gcalc_result_receiver *storage);
+  int get_line_result(res_point *cur, Gcalc_result_receiver *storage);
+
+  void free_result(res_point *res);
+};
+
+#endif /*GCALC_TOOLS_INCLUDED*/
+
diff --git a/sql/item_create.cc b/sql/item_create.cc
index b3a0e7cf3b2..8ff38bc4bc0 100644
--- a/sql/item_create.cc
+++ b/sql/item_create.cc
@@ -509,6 +509,19 @@ protected:
 
 
 #ifdef HAVE_SPATIAL
+class Create_func_mbr_contains : public Create_func_arg2
+{
+  public:
+    virtual Item *create_2_arg(THD *thd, Item *arg1, Item *arg2);
+
+    static Create_func_mbr_contains s_singleton;
+
+  protected:
+    Create_func_mbr_contains() {}
+    virtual ~Create_func_mbr_contains() {}
+};
+
+
 class Create_func_contains : public Create_func_arg2
 {
 public:
@@ -749,6 +762,19 @@ protected:
 
 
 #ifdef HAVE_SPATIAL
+class Create_func_mbr_disjoint : public Create_func_arg2
+{
+  public:
+    virtual Item *create_2_arg(THD *thd, Item *arg1, Item *arg2);
+
+    static Create_func_mbr_disjoint s_singleton;
+
+  protected:
+    Create_func_mbr_disjoint() {}
+    virtual ~Create_func_mbr_disjoint() {}
+};
+
+
 class Create_func_disjoint : public Create_func_arg2
 {
 public:
@@ -760,6 +786,19 @@ protected:
   Create_func_disjoint() {}
   virtual ~Create_func_disjoint() {}
 };
+
+
+class Create_func_distance : public Create_func_arg2
+{
+  public:
+    virtual Item* create_2_arg(THD *thd, Item *arg1, Item *arg2);
+
+    static Create_func_distance s_singleton;
+
+  protected:
+    Create_func_distance() {}
+    virtual ~Create_func_distance() {}
+};
 #endif
 
 
@@ -833,6 +872,19 @@ protected:
 
 
 #ifdef HAVE_SPATIAL
+class Create_func_mbr_equals : public Create_func_arg2
+{
+  public:
+    virtual Item *create_2_arg(THD *thd, Item *arg1, Item *arg2);
+
+    static Create_func_mbr_equals s_singleton;
+
+  protected:
+    Create_func_mbr_equals() {}
+    virtual ~Create_func_mbr_equals() {}
+};
+
+
 class Create_func_equals : public Create_func_arg2
 {
 public:
@@ -1161,6 +1213,19 @@ protected:
 
 
 #ifdef HAVE_SPATIAL
+class Create_func_mbr_intersects : public Create_func_arg2
+{
+  public:
+    virtual Item *create_2_arg(THD *thd, Item *arg1, Item *arg2);
+
+    static Create_func_mbr_intersects s_singleton;
+
+  protected:
+    Create_func_mbr_intersects() {}
+    virtual ~Create_func_mbr_intersects() {}
+};
+
+
 class Create_func_intersects : public Create_func_arg2
 {
 public:
@@ -1172,7 +1237,72 @@ protected:
   Create_func_intersects() {}
   virtual ~Create_func_intersects() {}
 };
-#endif
+
+
+class Create_func_intersection : public Create_func_arg2
+{
+public:
+  virtual Item* create_2_arg(THD *thd, Item *arg1, Item *arg2);
+
+  static Create_func_intersection s_singleton;
+
+protected:
+  Create_func_intersection() {}
+  virtual ~Create_func_intersection() {}
+};
+
+
+class Create_func_difference : public Create_func_arg2
+{
+public:
+  virtual Item* create_2_arg(THD *thd, Item *arg1, Item *arg2);
+
+  static Create_func_difference s_singleton;
+
+protected:
+  Create_func_difference() {}
+  virtual ~Create_func_difference() {}
+};
+
+
+class Create_func_union : public Create_func_arg2
+{
+public:
+  virtual Item* create_2_arg(THD *thd, Item *arg1, Item *arg2);
+
+  static Create_func_union s_singleton;
+
+protected:
+  Create_func_union() {}
+  virtual ~Create_func_union() {}
+};
+
+
+class Create_func_symdifference : public Create_func_arg2
+{
+public:
+  virtual Item* create_2_arg(THD *thd, Item *arg1, Item *arg2);
+
+  static Create_func_symdifference s_singleton;
+
+protected:
+  Create_func_symdifference() {}
+  virtual ~Create_func_symdifference() {}
+};
+
+
+class Create_func_buffer : public Create_func_arg2
+{
+public:
+  virtual Item* create_2_arg(THD *thd, Item *arg1, Item *arg2);
+
+  static Create_func_buffer s_singleton;
+
+protected:
+  Create_func_buffer() {}
+  virtual ~Create_func_buffer() {}
+};
+#endif /*HAVE_SPATIAL*/
 
 
 class Create_func_is_free_lock : public Create_func_arg1
@@ -1604,6 +1734,19 @@ protected:
 
 
 #ifdef HAVE_SPATIAL
+class Create_func_mbr_overlaps : public Create_func_arg2
+{
+  public:
+    virtual Item *create_2_arg(THD *thd, Item *arg1, Item *arg2);
+
+    static Create_func_mbr_overlaps s_singleton;
+
+  protected:
+    Create_func_mbr_overlaps() {}
+    virtual ~Create_func_mbr_overlaps() {}
+};
+
+
 class Create_func_overlaps : public Create_func_arg2
 {
 public:
@@ -2199,6 +2342,19 @@ protected:
 
 
 #ifdef HAVE_SPATIAL
+class Create_func_mbr_within : public Create_func_arg2
+{
+  public:
+    virtual Item *create_2_arg(THD *thd, Item *arg1, Item *arg2);
+
+    static Create_func_mbr_within s_singleton;
+
+  protected:
+    Create_func_mbr_within() {}
+    virtual ~Create_func_mbr_within() {}
+};
+
+
 class Create_func_within : public Create_func_arg2
 {
 public:
@@ -2890,6 +3046,16 @@ Create_func_connection_id::create_builder(THD *thd)
 
 
 #ifdef HAVE_SPATIAL
+Create_func_mbr_contains Create_func_mbr_contains::s_singleton;
+
+Item*
+Create_func_mbr_contains::create_2_arg(THD *thd, Item *arg1, Item *arg2)
+{
+  return new (thd->mem_root) Item_func_spatial_mbr_rel(arg1, arg2,
+      Item_func::SP_CONTAINS_FUNC);
+}
+
+
 Create_func_contains Create_func_contains::s_singleton;
 
 Item*
@@ -3122,6 +3288,16 @@ Create_func_dimension::create_1_arg(THD *thd, Item *arg1)
 
 
 #ifdef HAVE_SPATIAL
+Create_func_mbr_disjoint Create_func_mbr_disjoint::s_singleton;
+
+Item*
+Create_func_mbr_disjoint::create_2_arg(THD *thd, Item *arg1, Item *arg2)
+{
+  return new (thd->mem_root) Item_func_spatial_mbr_rel(arg1, arg2,
+      Item_func::SP_DISJOINT_FUNC);
+}
+
+
 Create_func_disjoint Create_func_disjoint::s_singleton;
 
 Item*
@@ -3130,6 +3306,15 @@ Create_func_disjoint::create_2_arg(THD *thd, Item *arg1, Item *arg2)
   return new (thd->mem_root) Item_func_spatial_rel(arg1, arg2,
                                                    Item_func::SP_DISJOINT_FUNC);
 }
+
+
+Create_func_distance Create_func_distance::s_singleton;
+
+Item*
+Create_func_distance::create_2_arg(THD *thd, Item *arg1, Item *arg2)
+{
+  return new (thd->mem_root) Item_func_distance(arg1, arg2);
+}
 #endif
 
 
@@ -3225,6 +3410,16 @@ Create_func_envelope::create_1_arg(THD *thd, Item *arg1)
 
 
 #ifdef HAVE_SPATIAL
+Create_func_mbr_equals Create_func_mbr_equals::s_singleton;
+
+Item*
+Create_func_mbr_equals::create_2_arg(THD *thd, Item *arg1, Item *arg2)
+{
+  return new (thd->mem_root) Item_func_spatial_mbr_rel(arg1, arg2,
+      Item_func::SP_EQUALS_FUNC);
+}
+
+
 Create_func_equals Create_func_equals::s_singleton;
 
 Item*
@@ -3620,6 +3815,16 @@ Create_func_interiorringn::create_2_arg(THD *thd, Item *arg1, Item *arg2)
 
 
 #ifdef HAVE_SPATIAL
+Create_func_mbr_intersects Create_func_mbr_intersects::s_singleton;
+
+Item*
+Create_func_mbr_intersects::create_2_arg(THD *thd, Item *arg1, Item *arg2)
+{
+  return new (thd->mem_root) Item_func_spatial_mbr_rel(arg1, arg2,
+      Item_func::SP_INTERSECTS_FUNC);
+}
+
+
 Create_func_intersects Create_func_intersects::s_singleton;
 
 Item*
@@ -3628,7 +3833,56 @@ Create_func_intersects::create_2_arg(THD *thd, Item *arg1, Item *arg2)
   return new (thd->mem_root) Item_func_spatial_rel(arg1, arg2,
                                                    Item_func::SP_INTERSECTS_FUNC);
 }
-#endif
+
+
+Create_func_intersection Create_func_intersection::s_singleton;
+
+Item*
+Create_func_intersection::create_2_arg(THD *thd, Item *arg1, Item *arg2)
+{
+  return new (thd->mem_root) Item_func_spatial_operation(arg1, arg2,
+                               Gcalc_function::op_intersection);
+}
+
+
+Create_func_difference Create_func_difference::s_singleton;
+
+Item*
+Create_func_difference::create_2_arg(THD *thd, Item *arg1, Item *arg2)
+{
+  return new (thd->mem_root) Item_func_spatial_operation(arg1, arg2,
+                               Gcalc_function::op_difference);
+}
+
+
+Create_func_union Create_func_union::s_singleton;
+
+Item*
+Create_func_union::create_2_arg(THD *thd, Item *arg1, Item *arg2)
+{
+  return new (thd->mem_root) Item_func_spatial_operation(arg1, arg2,
+                               Gcalc_function::op_union);
+}
+
+
+Create_func_symdifference Create_func_symdifference::s_singleton;
+
+Item*
+Create_func_symdifference::create_2_arg(THD *thd, Item *arg1, Item *arg2)
+{
+  return new (thd->mem_root) Item_func_spatial_operation(arg1, arg2,
+                               Gcalc_function::op_symdifference);
+}
+
+
+Create_func_buffer Create_func_buffer::s_singleton;
+
+Item*
+Create_func_buffer::create_2_arg(THD *thd, Item *arg1, Item *arg2)
+{
+  return new (thd->mem_root) Item_func_buffer(arg1, arg2);
+}
+#endif /*HAVE_SPATAI*/
 
 
 Create_func_is_free_lock Create_func_is_free_lock::s_singleton;
@@ -4088,6 +4342,16 @@ Create_func_ord::create_1_arg(THD *thd, Item *arg1)
 
 
 #ifdef HAVE_SPATIAL
+Create_func_mbr_overlaps Create_func_mbr_overlaps::s_singleton;
+
+Item*
+Create_func_mbr_overlaps::create_2_arg(THD *thd, Item *arg1, Item *arg2)
+{
+  return new (thd->mem_root) Item_func_spatial_mbr_rel(arg1, arg2,
+      Item_func::SP_OVERLAPS_FUNC);
+}
+
+
 Create_func_overlaps Create_func_overlaps::s_singleton;
 
 Item*
@@ -4640,6 +4904,16 @@ Create_func_weekofyear::create_1_arg(THD *thd, Item *arg1)
 
 
 #ifdef HAVE_SPATIAL
+Create_func_mbr_within Create_func_mbr_within::s_singleton;
+
+Item*
+Create_func_mbr_within::create_2_arg(THD *thd, Item *arg1, Item *arg2)
+{
+  return new (thd->mem_root) Item_func_spatial_mbr_rel(arg1, arg2,
+      Item_func::SP_WITHIN_FUNC);
+}
+
+
 Create_func_within Create_func_within::s_singleton;
 
 Item*
@@ -4773,6 +5047,7 @@ static Native_func_registry func_array[] =
   { { C_STRING_WITH_LEN("BIN") }, BUILDER(Create_func_bin)},
   { { C_STRING_WITH_LEN("BIT_COUNT") }, BUILDER(Create_func_bit_count)},
   { { C_STRING_WITH_LEN("BIT_LENGTH") }, BUILDER(Create_func_bit_length)},
+  { { C_STRING_WITH_LEN("BUFFER") }, GEOM_BUILDER(Create_func_buffer)},
   { { C_STRING_WITH_LEN("CEIL") }, BUILDER(Create_func_ceiling)},
   { { C_STRING_WITH_LEN("CEILING") }, BUILDER(Create_func_ceiling)},
   { { C_STRING_WITH_LEN("CENTROID") }, GEOM_BUILDER(Create_func_centroid)},
@@ -4800,13 +5075,13 @@ static Native_func_registry func_array[] =
   { { C_STRING_WITH_LEN("DES_DECRYPT") }, BUILDER(Create_func_des_decrypt)},
   { { C_STRING_WITH_LEN("DES_ENCRYPT") }, BUILDER(Create_func_des_encrypt)},
   { { C_STRING_WITH_LEN("DIMENSION") }, GEOM_BUILDER(Create_func_dimension)},
-  { { C_STRING_WITH_LEN("DISJOINT") }, GEOM_BUILDER(Create_func_disjoint)},
+  { { C_STRING_WITH_LEN("DISJOINT") }, GEOM_BUILDER(Create_func_mbr_disjoint)},
   { { C_STRING_WITH_LEN("ELT") }, BUILDER(Create_func_elt)},
   { { C_STRING_WITH_LEN("ENCODE") }, BUILDER(Create_func_encode)},
   { { C_STRING_WITH_LEN("ENCRYPT") }, BUILDER(Create_func_encrypt)},
   { { C_STRING_WITH_LEN("ENDPOINT") }, GEOM_BUILDER(Create_func_endpoint)},
   { { C_STRING_WITH_LEN("ENVELOPE") }, GEOM_BUILDER(Create_func_envelope)},
-  { { C_STRING_WITH_LEN("EQUALS") }, GEOM_BUILDER(Create_func_equals)},
+  { { C_STRING_WITH_LEN("EQUALS") }, GEOM_BUILDER(Create_func_mbr_equals)},
   { { C_STRING_WITH_LEN("EXP") }, BUILDER(Create_func_exp)},
   { { C_STRING_WITH_LEN("EXPORT_SET") }, BUILDER(Create_func_export_set)},
   { { C_STRING_WITH_LEN("EXTERIORRING") }, GEOM_BUILDER(Create_func_exteriorring)},
@@ -4837,7 +5112,7 @@ static Native_func_registry func_array[] =
   { { C_STRING_WITH_LEN("INET_NTOA") }, BUILDER(Create_func_inet_ntoa)},
   { { C_STRING_WITH_LEN("INSTR") }, BUILDER(Create_func_instr)},
   { { C_STRING_WITH_LEN("INTERIORRINGN") }, GEOM_BUILDER(Create_func_interiorringn)},
-  { { C_STRING_WITH_LEN("INTERSECTS") }, GEOM_BUILDER(Create_func_intersects)},
+  { { C_STRING_WITH_LEN("INTERSECTS") }, GEOM_BUILDER(Create_func_mbr_intersects)},
   { { C_STRING_WITH_LEN("ISCLOSED") }, GEOM_BUILDER(Create_func_isclosed)},
   { { C_STRING_WITH_LEN("ISEMPTY") }, GEOM_BUILDER(Create_func_isempty)},
   { { C_STRING_WITH_LEN("ISNULL") }, BUILDER(Create_func_isnull)},
@@ -4866,13 +5141,13 @@ static Native_func_registry func_array[] =
   { { C_STRING_WITH_LEN("MAKETIME") }, BUILDER(Create_func_maketime)},
   { { C_STRING_WITH_LEN("MAKE_SET") }, BUILDER(Create_func_make_set)},
   { { C_STRING_WITH_LEN("MASTER_POS_WAIT") }, BUILDER(Create_func_master_pos_wait)},
-  { { C_STRING_WITH_LEN("MBRCONTAINS") }, GEOM_BUILDER(Create_func_contains)},
-  { { C_STRING_WITH_LEN("MBRDISJOINT") }, GEOM_BUILDER(Create_func_disjoint)},
-  { { C_STRING_WITH_LEN("MBREQUAL") }, GEOM_BUILDER(Create_func_equals)},
-  { { C_STRING_WITH_LEN("MBRINTERSECTS") }, GEOM_BUILDER(Create_func_intersects)},
-  { { C_STRING_WITH_LEN("MBROVERLAPS") }, GEOM_BUILDER(Create_func_overlaps)},
+  { { C_STRING_WITH_LEN("MBRCONTAINS") }, GEOM_BUILDER(Create_func_mbr_contains)},
+  { { C_STRING_WITH_LEN("MBRDISJOINT") }, GEOM_BUILDER(Create_func_mbr_disjoint)},
+  { { C_STRING_WITH_LEN("MBREQUAL") }, GEOM_BUILDER(Create_func_mbr_equals)},
+  { { C_STRING_WITH_LEN("MBRINTERSECTS") }, GEOM_BUILDER(Create_func_mbr_intersects)},
+  { { C_STRING_WITH_LEN("MBROVERLAPS") }, GEOM_BUILDER(Create_func_mbr_overlaps)},
   { { C_STRING_WITH_LEN("MBRTOUCHES") }, GEOM_BUILDER(Create_func_touches)},
-  { { C_STRING_WITH_LEN("MBRWITHIN") }, GEOM_BUILDER(Create_func_within)},
+  { { C_STRING_WITH_LEN("MBRWITHIN") }, GEOM_BUILDER(Create_func_mbr_within)},
   { { C_STRING_WITH_LEN("MD5") }, BUILDER(Create_func_md5)},
   { { C_STRING_WITH_LEN("MLINEFROMTEXT") }, GEOM_BUILDER(Create_func_geometry_from_text)},
   { { C_STRING_WITH_LEN("MLINEFROMWKB") }, GEOM_BUILDER(Create_func_geometry_from_wkb)},
@@ -4895,7 +5170,7 @@ static Native_func_registry func_array[] =
   { { C_STRING_WITH_LEN("OCT") }, BUILDER(Create_func_oct)},
   { { C_STRING_WITH_LEN("OCTET_LENGTH") }, BUILDER(Create_func_length)},
   { { C_STRING_WITH_LEN("ORD") }, BUILDER(Create_func_ord)},
-  { { C_STRING_WITH_LEN("OVERLAPS") }, GEOM_BUILDER(Create_func_overlaps)},
+  { { C_STRING_WITH_LEN("OVERLAPS") }, GEOM_BUILDER(Create_func_mbr_overlaps)},
   { { C_STRING_WITH_LEN("PERIOD_ADD") }, BUILDER(Create_func_period_add)},
   { { C_STRING_WITH_LEN("PERIOD_DIFF") }, BUILDER(Create_func_period_diff)},
   { { C_STRING_WITH_LEN("PI") }, BUILDER(Create_func_pi)},
@@ -4930,6 +5205,64 @@ static Native_func_registry func_array[] =
   { { C_STRING_WITH_LEN("STARTPOINT") }, GEOM_BUILDER(Create_func_startpoint)},
   { { C_STRING_WITH_LEN("STRCMP") }, BUILDER(Create_func_strcmp)},
   { { C_STRING_WITH_LEN("STR_TO_DATE") }, BUILDER(Create_func_str_to_date)},
+  { { C_STRING_WITH_LEN("ST_AREA") }, GEOM_BUILDER(Create_func_area)},
+  { { C_STRING_WITH_LEN("ST_ASBINARY") }, GEOM_BUILDER(Create_func_as_wkb)},
+  { { C_STRING_WITH_LEN("ST_ASTEXT") }, GEOM_BUILDER(Create_func_as_wkt)},
+  { { C_STRING_WITH_LEN("ST_ASWKB") }, GEOM_BUILDER(Create_func_as_wkb)},
+  { { C_STRING_WITH_LEN("ST_ASWKT") }, GEOM_BUILDER(Create_func_as_wkt)},
+  { { C_STRING_WITH_LEN("ST_BUFFER") }, GEOM_BUILDER(Create_func_buffer)},
+  { { C_STRING_WITH_LEN("ST_CENTROID") }, GEOM_BUILDER(Create_func_centroid)},
+  { { C_STRING_WITH_LEN("ST_CONTAINS") }, GEOM_BUILDER(Create_func_contains)},
+  { { C_STRING_WITH_LEN("ST_CROSSES") }, GEOM_BUILDER(Create_func_crosses)},
+  { { C_STRING_WITH_LEN("ST_DIFFERENCE") }, GEOM_BUILDER(Create_func_difference)},
+  { { C_STRING_WITH_LEN("ST_DIMENSION") }, GEOM_BUILDER(Create_func_dimension)},
+  { { C_STRING_WITH_LEN("ST_DISJOINT") }, GEOM_BUILDER(Create_func_disjoint)},
+  { { C_STRING_WITH_LEN("ST_DISTANCE") }, GEOM_BUILDER(Create_func_distance)},
+  { { C_STRING_WITH_LEN("ST_ENDPOINT") }, GEOM_BUILDER(Create_func_endpoint)},
+  { { C_STRING_WITH_LEN("ST_ENVELOPE") }, GEOM_BUILDER(Create_func_envelope)},
+  { { C_STRING_WITH_LEN("ST_EQUALS") }, GEOM_BUILDER(Create_func_mbr_equals)},
+  { { C_STRING_WITH_LEN("ST_EXTERIORRING") }, GEOM_BUILDER(Create_func_exteriorring)},
+  { { C_STRING_WITH_LEN("ST_GEOMCOLLFROMTEXT") }, GEOM_BUILDER(Create_func_geometry_from_text)},
+  { { C_STRING_WITH_LEN("ST_GEOMCOLLFROMWKB") }, GEOM_BUILDER(Create_func_geometry_from_wkb)},
+  { { C_STRING_WITH_LEN("ST_GEOMETRYCOLLECTIONFROMTEXT") }, GEOM_BUILDER(Create_func_geometry_from_text)},
+  { { C_STRING_WITH_LEN("ST_GEOMETRYCOLLECTIONFROMWKB") }, GEOM_BUILDER(Create_func_geometry_from_wkb)},
+  { { C_STRING_WITH_LEN("ST_GEOMETRYFROMTEXT") }, GEOM_BUILDER(Create_func_geometry_from_text)},
+  { { C_STRING_WITH_LEN("ST_GEOMETRYFROMWKB") }, GEOM_BUILDER(Create_func_geometry_from_wkb)},
+  { { C_STRING_WITH_LEN("ST_GEOMETRYN") }, GEOM_BUILDER(Create_func_geometryn)},
+  { { C_STRING_WITH_LEN("ST_GEOMETRYTYPE") }, GEOM_BUILDER(Create_func_geometry_type)},
+  { { C_STRING_WITH_LEN("ST_GEOMFROMTEXT") }, GEOM_BUILDER(Create_func_geometry_from_text)},
+  { { C_STRING_WITH_LEN("ST_GEOMFROMWKB") }, GEOM_BUILDER(Create_func_geometry_from_wkb)},
+  { { C_STRING_WITH_LEN("ST_EQUALS") }, GEOM_BUILDER(Create_func_equals)},
+  { { C_STRING_WITH_LEN("ST_INTERIORRINGN") }, GEOM_BUILDER(Create_func_interiorringn)},
+  { { C_STRING_WITH_LEN("ST_INTERSECTS") }, GEOM_BUILDER(Create_func_intersects)},
+  { { C_STRING_WITH_LEN("ST_INTERSECTION") }, GEOM_BUILDER(Create_func_intersection)},
+  { { C_STRING_WITH_LEN("ST_ISCLOSED") }, GEOM_BUILDER(Create_func_isclosed)},
+  { { C_STRING_WITH_LEN("ST_ISEMPTY") }, GEOM_BUILDER(Create_func_isempty)},
+  { { C_STRING_WITH_LEN("ST_ISSIMPLE") }, GEOM_BUILDER(Create_func_issimple)},
+  { { C_STRING_WITH_LEN("ST_LENGTH") }, GEOM_BUILDER(Create_func_glength)},
+  { { C_STRING_WITH_LEN("ST_LINEFROMTEXT") }, GEOM_BUILDER(Create_func_geometry_from_text)},
+  { { C_STRING_WITH_LEN("ST_LINEFROMWKB") }, GEOM_BUILDER(Create_func_geometry_from_wkb)},
+  { { C_STRING_WITH_LEN("ST_LINESTRINGFROMTEXT") }, GEOM_BUILDER(Create_func_geometry_from_text)},
+  { { C_STRING_WITH_LEN("ST_LINESTRINGFROMWKB") }, GEOM_BUILDER(Create_func_geometry_from_wkb)},
+  { { C_STRING_WITH_LEN("ST_NUMGEOMETRIES") }, GEOM_BUILDER(Create_func_numgeometries)},
+  { { C_STRING_WITH_LEN("ST_NUMINTERIORRINGS") }, GEOM_BUILDER(Create_func_numinteriorring)},
+  { { C_STRING_WITH_LEN("ST_NUMPOINTS") }, GEOM_BUILDER(Create_func_numpoints)},
+  { { C_STRING_WITH_LEN("ST_OVERLAPS") }, GEOM_BUILDER(Create_func_overlaps)},
+  { { C_STRING_WITH_LEN("ST_POINTFROMTEXT") }, GEOM_BUILDER(Create_func_geometry_from_text)},
+  { { C_STRING_WITH_LEN("ST_POINTFROMWKB") }, GEOM_BUILDER(Create_func_geometry_from_wkb)},
+  { { C_STRING_WITH_LEN("ST_POINTN") }, GEOM_BUILDER(Create_func_pointn)},
+  { { C_STRING_WITH_LEN("ST_POLYFROMTEXT") }, GEOM_BUILDER(Create_func_geometry_from_text)},
+  { { C_STRING_WITH_LEN("ST_POLYFROMWKB") }, GEOM_BUILDER(Create_func_geometry_from_wkb)},
+  { { C_STRING_WITH_LEN("ST_POLYGONFROMTEXT") }, GEOM_BUILDER(Create_func_geometry_from_text)},
+  { { C_STRING_WITH_LEN("ST_POLYGONFROMWKB") }, GEOM_BUILDER(Create_func_geometry_from_wkb)},
+  { { C_STRING_WITH_LEN("ST_SRID") }, GEOM_BUILDER(Create_func_srid)},
+  { { C_STRING_WITH_LEN("ST_STARTPOINT") }, GEOM_BUILDER(Create_func_startpoint)},
+  { { C_STRING_WITH_LEN("ST_SYMDIFFERENCE") }, GEOM_BUILDER(Create_func_symdifference)},
+  { { C_STRING_WITH_LEN("ST_TOUCHES") }, GEOM_BUILDER(Create_func_touches)},
+  { { C_STRING_WITH_LEN("ST_UNION") }, GEOM_BUILDER(Create_func_union)},
+  { { C_STRING_WITH_LEN("ST_WITHIN") }, GEOM_BUILDER(Create_func_mbr_within)},
+  { { C_STRING_WITH_LEN("ST_X") }, GEOM_BUILDER(Create_func_x)},
+  { { C_STRING_WITH_LEN("ST_Y") }, GEOM_BUILDER(Create_func_y)},
   { { C_STRING_WITH_LEN("SUBSTRING_INDEX") }, BUILDER(Create_func_substr_index)},
   { { C_STRING_WITH_LEN("SUBTIME") }, BUILDER(Create_func_subtime)},
   { { C_STRING_WITH_LEN("TAN") }, BUILDER(Create_func_tan)},
diff --git a/sql/item_geofunc.cc b/sql/item_geofunc.cc
index 8c38cb2a859..9432de95182 100644
--- a/sql/item_geofunc.cc
+++ b/sql/item_geofunc.cc
@@ -1,4 +1,4 @@
-/* Copyright (C) 2003-2006 MySQL AB
+/* Copyright (c) 2003, 2010, Oracle and/or its affiliates. All rights reserved.
 
    This program is free software; you can redistribute it and/or modify
    it under the terms of the GNU General Public License as published by
@@ -10,8 +10,8 @@
    GNU General Public License for more details.
 
    You should have received a copy of the GNU General Public License
-   along with this program; if not, write to the Free Software
-   Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA */
+   along with this program; if not, write to the Free Software Foundation,
+   51 Franklin Street, Suite 500, Boston, MA 02110-1335 USA */
 
 
 /**
@@ -26,9 +26,15 @@
 #endif
 
 #include "mysql_priv.h"
+/*
+  It is necessary to include set_var.h instead of item.h because there
+  are dependencies on include order for set_var.h and item.h. This
+  will be resolved later.
+*/
 #ifdef HAVE_SPATIAL
 #include <m_ctype.h>
 
+
 Field *Item_geometry_func::tmp_table_field(TABLE *t_arg)
 {
   Field *result;
@@ -42,7 +48,7 @@ void Item_geometry_func::fix_length_and_dec()
 {
   collation.set(&my_charset_bin);
   decimals=0;
-  max_length= max_field_size;
+  max_length= (uint32) 4294967295U;
   maybe_null= 1;
 }
 
@@ -134,6 +140,7 @@ String *Item_func_as_wkt::val_str(String *str)
 
 void Item_func_as_wkt::fix_length_and_dec()
 {
+  collation.set(default_charset(), DERIVATION_COERCIBLE, MY_REPERTOIRE_ASCII);
   max_length=MAX_BLOB_WIDTH;
   maybe_null= 1;
 }
@@ -360,8 +367,8 @@ String *Item_func_point::val_str(String *str)
   uint32 srid= 0;
 
   if ((null_value= (args[0]->null_value ||
-		    args[1]->null_value ||
-                    str->realloc(4/*SRID*/ + 1 + 4 + SIZEOF_STORED_DOUBLE*2))))
+                    args[1]->null_value ||
+                    str->realloc(4/*SRID*/ + 1 + 4 + SIZEOF_STORED_DOUBLE * 2))))
     return 0;
 
   str->set_charset(&my_charset_bin);
@@ -508,7 +515,33 @@ err:
   Functions for spatial relations
 */
 
-longlong Item_func_spatial_rel::val_int()
+const char *Item_func_spatial_mbr_rel::func_name() const 
+{ 
+  switch (spatial_rel) {
+    case SP_CONTAINS_FUNC:
+      return "mbrcontains";
+    case SP_WITHIN_FUNC:
+      return "mbrwithin";
+    case SP_EQUALS_FUNC:
+      return "mbrequals";
+    case SP_DISJOINT_FUNC:
+      return "mbrdisjoint";
+    case SP_INTERSECTS_FUNC:
+      return "mbrintersects";
+    case SP_TOUCHES_FUNC:
+      return "mbrtouches";
+    case SP_CROSSES_FUNC:
+      return "mbrcrosses";
+    case SP_OVERLAPS_FUNC:
+      return "mbroverlaps";
+    default:
+      DBUG_ASSERT(0);  // Should never happened
+      return "mbrsp_unknown"; 
+  }
+}
+
+
+longlong Item_func_spatial_mbr_rel::val_int()
 {
   DBUG_ASSERT(fixed == 1);
   String *res1= args[0]->val_str(&cmp.value1);
@@ -553,6 +586,867 @@ longlong Item_func_spatial_rel::val_int()
 }
 
 
+Item_func_spatial_rel::Item_func_spatial_rel(Item *a,Item *b,
+                                             enum Functype sp_rel) :
+    Item_bool_func2(a,b), collector()
+{
+  spatial_rel = sp_rel;
+}
+
+
+Item_func_spatial_rel::~Item_func_spatial_rel()
+{
+}
+
+
+const char *Item_func_spatial_rel::func_name() const 
+{ 
+  switch (spatial_rel) {
+    case SP_CONTAINS_FUNC:
+      return "st_contains";
+    case SP_WITHIN_FUNC:
+      return "st_within";
+    case SP_EQUALS_FUNC:
+      return "st_equals";
+    case SP_DISJOINT_FUNC:
+      return "st_disjoint";
+    case SP_INTERSECTS_FUNC:
+      return "st_intersects";
+    case SP_TOUCHES_FUNC:
+      return "st_touches";
+    case SP_CROSSES_FUNC:
+      return "st_crosses";
+    case SP_OVERLAPS_FUNC:
+      return "st_overlaps";
+    default:
+      DBUG_ASSERT(0);  // Should never happened
+      return "sp_unknown"; 
+  }
+}
+
+
+static double count_edge_t(const Gcalc_heap::Info *ea,
+                           const Gcalc_heap::Info *eb,
+                           const Gcalc_heap::Info *v,
+                           double &ex, double &ey, double &vx, double &vy,
+                           double &e_sqrlen)
+{
+  ex= eb->x - ea->x;
+  ey= eb->y - ea->y;
+  vx= v->x - ea->x;
+  vy= v->y - ea->y;
+  e_sqrlen= ex * ex + ey * ey;
+  return (ex * vx + ey * vy) / e_sqrlen;
+}
+
+
+static double distance_to_line(double ex, double ey, double vx, double vy,
+                               double e_sqrlen)
+{
+  return fabs(vx * ey - vy * ex) / sqrt(e_sqrlen);
+}
+
+
+static double distance_points(const Gcalc_heap::Info *a,
+                              const Gcalc_heap::Info *b)
+{
+  double x= a->x - b->x;
+  double y= a->y - b->y;
+  return sqrt(x * x + y * y);
+}
+
+
+/*
+  Calculates the distance between objects.
+*/
+
+static int calc_distance(double *result, Gcalc_heap *collector, uint obj2_si,
+                         Gcalc_function *func, Gcalc_scan_iterator *scan_it)
+{
+  bool above_cur_point, cur_point_edge;
+  const Gcalc_scan_iterator::point *evpos;
+  const Gcalc_heap::Info *cur_point, *dist_point;
+  Gcalc_scan_events ev;
+  double t, distance, cur_distance;
+  double ex, ey, vx, vy, e_sqrlen;
+
+  DBUG_ENTER("calc_distance");
+
+  above_cur_point= false;
+  distance= DBL_MAX;
+
+  while (scan_it->more_points())
+  {
+    if (scan_it->step())
+      goto mem_error;
+    evpos= scan_it->get_event_position();
+    ev= scan_it->get_event();
+    cur_point= evpos->pi;
+
+    /*
+       handling intersection we only need to check if it's the intersecion
+       of objects 1 and 2. In this case distance is 0
+    */
+    if (ev == scev_intersection)
+    {
+      if ((evpos->get_next()->pi->shape >= obj2_si) !=
+            (cur_point->shape >= obj2_si))
+      {
+        distance= 0;
+        goto exit;
+      }
+      continue;
+    }
+
+    /*
+       if we get 'scev_point | scev_end | scev_two_ends' we don't need
+       to check for intersection of objects.
+       Though we need to calculate distances.
+    */
+    if (ev & (scev_point | scev_end | scev_two_ends))
+      goto calculate_distance;
+
+    goto calculate_distance;
+    /*
+       having these events we need to check for possible intersection
+       of objects
+       scev_thread | scev_two_threads | scev_single_point
+    */
+    DBUG_ASSERT(ev & (scev_thread | scev_two_threads | scev_single_point));
+
+    func->clear_state();
+    for (Gcalc_point_iterator pit(scan_it); pit.point() != evpos; ++pit)
+    {
+      gcalc_shape_info si= pit.point()->get_shape();
+      if ((func->get_shape_kind(si) == Gcalc_function::shape_polygon))
+        func->invert_state(si);
+    }
+    func->invert_state(evpos->get_shape());
+    if (func->count())
+    {
+      /* Point of one object is inside the other - intersection found */
+      distance= 0;
+      goto exit;
+    }
+
+
+calculate_distance:
+    if (cur_point->shape >= obj2_si)
+      continue;
+    cur_point_edge= !cur_point->is_bottom();
+
+    for (dist_point= collector->get_first(); dist_point; dist_point= dist_point->get_next())
+    {
+      /* We only check vertices of object 2 */
+      if (dist_point->shape < obj2_si)
+        continue;
+
+      /* if we have an edge to check */
+      if (dist_point->left)
+      {
+        t= count_edge_t(dist_point, dist_point->left, cur_point,
+                        ex, ey, vx, vy, e_sqrlen);
+        if ((t > 0.0) && (t < 1.0))
+        {
+          cur_distance= distance_to_line(ex, ey, vx, vy, e_sqrlen);
+          if (distance > cur_distance)
+            distance= cur_distance;
+        }
+      }
+      if (cur_point_edge)
+      {
+        t= count_edge_t(cur_point, cur_point->left, dist_point,
+                        ex, ey, vx, vy, e_sqrlen);
+        if ((t > 0.0) && (t < 1.0))
+        {
+          cur_distance= distance_to_line(ex, ey, vx, vy, e_sqrlen);
+          if (distance > cur_distance)
+            distance= cur_distance;
+        }
+      }
+      cur_distance= distance_points(cur_point, dist_point);
+      if (distance > cur_distance)
+        distance= cur_distance;
+    }
+  }
+
+exit:
+  *result= distance;
+  DBUG_RETURN(0);
+
+mem_error:
+  DBUG_RETURN(1);
+}
+
+
+#define GIS_ZERO 0.00000000001
+
+int Item_func_spatial_rel::func_touches()
+{
+  bool above_cur_point;
+  double x1, x2, y1, y2, ex, ey;
+  double distance, area;
+  int result= 0;
+  int cur_func= 0;
+
+  Gcalc_operation_transporter trn(&func, &collector);
+
+  String *res1= args[0]->val_str(&tmp_value1);
+  String *res2= args[1]->val_str(&tmp_value2);
+  Geometry_buffer buffer1, buffer2;
+  Geometry *g1, *g2;
+  int obj2_si;
+
+  DBUG_ENTER("Item_func_spatial_rel::func_touches");
+  DBUG_ASSERT(fixed == 1);
+
+  if ((null_value= (args[0]->null_value || args[1]->null_value ||
+          !(g1= Geometry::construct(&buffer1, res1->ptr(), res1->length())) ||
+          !(g2= Geometry::construct(&buffer2, res2->ptr(), res2->length())))))
+    goto mem_error;
+
+  if ((g1->get_class_info()->m_type_id == Geometry::wkb_point) &&
+      (g2->get_class_info()->m_type_id == Geometry::wkb_point))
+  {
+    if (((Gis_point *) g1)->get_xy(&x1, &y1) ||
+        ((Gis_point *) g2)->get_xy(&x2, &y2))
+      goto mem_error;
+    ex= x2 - x1;
+    ey= y2 - y1;
+    DBUG_RETURN((ex * ex + ey * ey) < GIS_ZERO);
+  }
+
+  if (func.reserve_op_buffer(1))
+    goto mem_error;
+  func.add_operation(Gcalc_function::op_intersection, 2);
+
+  if (g1->store_shapes(&trn))
+    goto mem_error;
+  obj2_si= func.get_nshapes();
+
+  if (g2->store_shapes(&trn) || func.alloc_states())
+    goto mem_error;
+
+  collector.prepare_operation();
+  scan_it.init(&collector);
+
+  if (calc_distance(&distance, &collector, obj2_si, &func, &scan_it))
+   goto mem_error;
+  if (distance > GIS_ZERO)
+    goto exit;
+
+  scan_it.reset();
+  scan_it.init(&collector);
+
+  above_cur_point= false;
+  distance= DBL_MAX;
+
+  while (scan_it.more_trapezoids())
+  {
+    if (scan_it.step())
+      goto mem_error;
+
+    func.clear_state();
+    for (Gcalc_trapezoid_iterator ti(&scan_it); ti.more(); ++ti)
+    {
+      gcalc_shape_info si= ti.lb()->get_shape();
+      if ((func.get_shape_kind(si) == Gcalc_function::shape_polygon))
+      {
+        func.invert_state(si);
+        cur_func= func.count();
+      }
+      if (cur_func)
+      {
+        area= scan_it.get_h() *
+              ((ti.rb()->x - ti.lb()->x) + (ti.rt()->x - ti.lt()->x));
+        if (area > GIS_ZERO)
+        {
+          result= 0;
+          goto exit;
+        }
+      }
+    }
+  }
+  result= 1;
+
+exit:
+  collector.reset();
+  func.reset();
+  scan_it.reset();
+  DBUG_RETURN(result);
+mem_error:
+  null_value= 1;
+  DBUG_RETURN(0);
+}
+
+
+int Item_func_spatial_rel::func_equals()
+{
+  Gcalc_heap::Info *pi_s1, *pi_s2;
+  Gcalc_heap::Info *cur_pi= collector.get_first();
+  double d;
+
+  if (!cur_pi)
+    return 1;
+
+  do {
+    pi_s1= cur_pi;
+    pi_s2= 0;
+    while ((cur_pi= cur_pi->get_next()))
+    {
+      d= fabs(pi_s1->x - cur_pi->x) + fabs(pi_s1->y - cur_pi->y);
+      if (d > GIS_ZERO)
+        break;
+      if (!pi_s2 && pi_s1->shape != cur_pi->shape)
+        pi_s2= cur_pi;
+    }
+
+    if (!pi_s2)
+      return 0;
+  } while (cur_pi);
+
+  return 1;
+}
+
+
+longlong Item_func_spatial_rel::val_int()
+{
+  DBUG_ENTER("Item_func_spatial_rel::val_int");
+  DBUG_ASSERT(fixed == 1);
+  String *res1;
+  String *res2;
+  Geometry_buffer buffer1, buffer2;
+  Geometry *g1, *g2;
+  int result= 0;
+  int mask= 0;
+
+  if (spatial_rel == SP_TOUCHES_FUNC)
+    DBUG_RETURN(func_touches());
+
+  res1= args[0]->val_str(&tmp_value1);
+  res2= args[1]->val_str(&tmp_value2);
+  Gcalc_operation_transporter trn(&func, &collector);
+
+  if (func.reserve_op_buffer(1))
+    DBUG_RETURN(0);
+
+  switch (spatial_rel) {
+    case SP_CONTAINS_FUNC:
+      mask= 1;
+      func.add_operation(Gcalc_function::op_backdifference, 2);
+      break;
+    case SP_WITHIN_FUNC:
+      mask= 1;
+      func.add_operation(Gcalc_function::op_difference, 2);
+      break;
+    case SP_EQUALS_FUNC:
+      break;
+    case SP_DISJOINT_FUNC:
+      mask= 1;
+      func.add_operation(Gcalc_function::op_intersection, 2);
+      break;
+    case SP_INTERSECTS_FUNC:
+      func.add_operation(Gcalc_function::op_intersection, 2);
+      break;
+    case SP_OVERLAPS_FUNC:
+      func.add_operation(Gcalc_function::op_backdifference, 2);
+      break;
+    case SP_CROSSES_FUNC:
+      func.add_operation(Gcalc_function::op_intersection, 2);
+      break;
+    default:
+      DBUG_ASSERT(FALSE);
+      break;
+  }
+
+
+  if ((null_value=
+       (args[0]->null_value || args[1]->null_value ||
+	!(g1= Geometry::construct(&buffer1, res1->ptr(), res1->length())) ||
+	!(g2= Geometry::construct(&buffer2, res2->ptr(), res2->length())) ||
+	g1->store_shapes(&trn) || g2->store_shapes(&trn))))
+    goto exit;
+
+  collector.prepare_operation();
+  scan_it.init(&collector);
+  if (spatial_rel == SP_EQUALS_FUNC)
+  {
+    result= (g1->get_class_info()->m_type_id == g1->get_class_info()->m_type_id) &&
+            func_equals();
+    goto exit;
+  }
+
+  if (func.alloc_states())
+    goto exit;
+
+  result= func.find_function(scan_it) ^ mask;
+
+exit:
+  collector.reset();
+  func.reset();
+  scan_it.reset();
+  DBUG_RETURN(result);
+}
+
+
+Item_func_spatial_operation::~Item_func_spatial_operation()
+{
+}
+
+
+String *Item_func_spatial_operation::val_str(String *str_value)
+{
+  DBUG_ENTER("Item_func_spatial_operation::val_str");
+  DBUG_ASSERT(fixed == 1);
+  String *res1= args[0]->val_str(&tmp_value1);
+  String *res2= args[1]->val_str(&tmp_value2);
+  Geometry_buffer buffer1, buffer2;
+  Geometry *g1, *g2;
+  uint32 srid= 0;
+  Gcalc_operation_transporter trn(&func, &collector);
+
+  if (func.reserve_op_buffer(1))
+    DBUG_RETURN(0);
+  func.add_operation(spatial_op, 2);
+
+  if ((null_value=
+       (args[0]->null_value || args[1]->null_value ||
+	!(g1= Geometry::construct(&buffer1, res1->ptr(), res1->length())) ||
+	!(g2= Geometry::construct(&buffer2, res2->ptr(), res2->length())) ||
+	g1->store_shapes(&trn) || g2->store_shapes(&trn))))
+    goto exit;
+
+  
+  collector.prepare_operation();
+  scan_it.init(&collector);
+  if (func.alloc_states())
+    goto exit;
+
+  operation.init(&func);
+
+  if (operation.count_all(&collector) ||
+      operation.get_result(&res_receiver))
+    goto exit;
+
+
+  str_value->set_charset(&my_charset_bin);
+  if (str_value->reserve(SRID_SIZE, 512))
+    goto exit;
+  str_value->length(0);
+  str_value->q_append(srid);
+
+  if (!Geometry::create_from_opresult(&buffer1, str_value, res_receiver))
+    goto exit;
+
+exit:
+  collector.reset();
+  func.reset();
+  scan_it.reset();
+  res_receiver.reset();
+  DBUG_RETURN(str_value);
+}
+
+
+const char *Item_func_spatial_operation::func_name() const
+{ 
+  switch (spatial_op) {
+    case Gcalc_function::op_intersection:
+      return "st_intersection";
+    case Gcalc_function::op_difference:
+      return "st_difference";
+    case Gcalc_function::op_union:
+      return "st_union";
+    case Gcalc_function::op_symdifference:
+      return "st_symdifference";
+    default:
+      DBUG_ASSERT(0);  // Should never happen
+      return "sp_unknown"; 
+  }
+}
+
+
+static const int SINUSES_CALCULATED= 32;
+static double n_sinus[SINUSES_CALCULATED+1]=
+{
+  0,
+  0.04906767432741802,
+  0.0980171403295606,
+  0.1467304744553618,
+  0.1950903220161283,
+  0.2429801799032639,
+  0.2902846772544623,
+  0.3368898533922201,
+  0.3826834323650898,
+  0.4275550934302821,
+  0.4713967368259976,
+  0.5141027441932217,
+  0.5555702330196022,
+  0.5956993044924334,
+  0.6343932841636455,
+  0.6715589548470183,
+  0.7071067811865475,
+  0.7409511253549591,
+  0.773010453362737,
+  0.8032075314806448,
+  0.8314696123025452,
+  0.8577286100002721,
+  0.8819212643483549,
+  0.9039892931234433,
+  0.9238795325112867,
+  0.9415440651830208,
+  0.9569403357322089,
+  0.970031253194544,
+  0.9807852804032304,
+  0.989176509964781,
+  0.9951847266721968,
+  0.9987954562051724,
+  1
+};
+
+
+static void get_n_sincos(int n, double *sinus, double *cosinus)
+{
+  DBUG_ASSERT(n > 0 && n < SINUSES_CALCULATED*2+1);
+  if (n < (SINUSES_CALCULATED + 1))
+  {
+    *sinus= n_sinus[n];
+    *cosinus= n_sinus[SINUSES_CALCULATED - n];
+  }
+  else
+  {
+    n-= SINUSES_CALCULATED;
+    *sinus= n_sinus[SINUSES_CALCULATED - n];
+    *cosinus= -n_sinus[n];
+  }
+}
+
+
+static int fill_half_circle(Gcalc_shape_transporter *trn, double x, double y,
+                            double ax, double ay)
+{
+  double n_sin, n_cos;
+  double x_n, y_n;
+  for (int n = 1; n < (SINUSES_CALCULATED * 2 - 1); n++)
+  {
+    get_n_sincos(n, &n_sin, &n_cos);
+    x_n= ax * n_cos - ay * n_sin;
+    y_n= ax * n_sin + ay * n_cos;
+    if (trn->add_point(x_n + x, y_n + y))
+      return 1;
+  }
+  return 0;
+}
+
+
+static int fill_gap(Gcalc_shape_transporter *trn,
+                    double x, double y,
+                    double ax, double ay, double bx, double by, double d,
+                    bool *empty_gap)
+{
+  double ab= ax * bx + ay * by;
+  double cosab= ab / (d * d) + GIS_ZERO;
+  double n_sin, n_cos;
+  double x_n, y_n;
+  int n=1;
+
+  *empty_gap= true;
+  for (;;)
+  {
+    get_n_sincos(n++, &n_sin, &n_cos);
+    if (n_cos <= cosab)
+      break;
+    *empty_gap= false;
+    x_n= ax * n_cos - ay * n_sin;
+    y_n= ax * n_sin + ay * n_cos;
+    if (trn->add_point(x_n + x, y_n + y))
+      return 1;
+  }
+  return 0;
+}
+
+
+/*
+  Calculates the vector (p2,p1) and
+  negatively orthogonal to it with the length of d.
+  The result is (ex,ey) - the vector, (px,py) - the orthogonal.
+*/
+
+static void calculate_perpendicular(
+    double x1, double y1, double x2, double y2, double d,
+    double *ex, double *ey,
+    double *px, double *py)
+{
+  double q;
+  *ex= x1 - x2;
+  *ey= y1 - y2;
+  q= d / sqrt((*ex) * (*ex) + (*ey) * (*ey));
+  *px= (*ey) * q;
+  *py= -(*ex) * q;
+}
+
+
+int Item_func_buffer::Transporter::single_point(double x, double y)
+{
+  return add_point_buffer(x, y);
+}
+
+
+int Item_func_buffer::Transporter::add_edge_buffer(
+  double x3, double y3, bool round_p1, bool round_p2)
+{
+  Gcalc_operation_transporter trn(m_fn, m_heap);
+  double e1_x, e1_y, e2_x, e2_y, p1_x, p1_y, p2_x, p2_y;
+  double e1e2;
+  double sin1, cos1;
+  double x_n, y_n;
+  bool empty_gap1, empty_gap2;
+
+  ++m_nshapes;
+  if (trn.start_simple_poly())
+    return 1;
+
+  calculate_perpendicular(x1, y1, x2, y2, m_d, &e1_x, &e1_y, &p1_x, &p1_y);
+  calculate_perpendicular(x3, y3, x2, y2, m_d, &e2_x, &e2_y, &p2_x, &p2_y);
+
+  e1e2= e1_x * e2_y - e2_x * e1_y;
+  sin1= n_sinus[1];
+  cos1= n_sinus[31];
+  if (e1e2 < 0)
+  {
+    empty_gap2= false;
+    x_n= x2 + p2_x * cos1 - p2_y * sin1;
+    y_n= y2 + p2_y * cos1 + p2_x * sin1;
+    if (fill_gap(&trn, x2, y2, -p1_x,-p1_y, p2_x,p2_y, m_d, &empty_gap1) ||
+        trn.add_point(x2 + p2_x, y2 + p2_y) ||
+        trn.add_point(x_n, y_n))
+      return 1;
+  }
+  else
+  {
+    x_n= x2 - p2_x * cos1 - p2_y * sin1;
+    y_n= y2 - p2_y * cos1 + p2_x * sin1;
+    if (trn.add_point(x_n, y_n) ||
+        trn.add_point(x2 - p2_x, y2 - p2_y) ||
+        fill_gap(&trn, x2, y2, -p2_x, -p2_y, p1_x, p1_y, m_d, &empty_gap2))
+      return 1;
+    empty_gap1= false;
+  }
+  if ((!empty_gap2 && trn.add_point(x2 + p1_x, y2 + p1_y)) ||
+      trn.add_point(x1 + p1_x, y1 + p1_y))
+    return 1;
+
+  if (round_p1 && fill_half_circle(&trn, x1, y1, p1_x, p1_y))
+    return 1;
+
+  if (trn.add_point(x1 - p1_x, y1 - p1_y) ||
+      (!empty_gap1 && trn.add_point(x2 - p1_x, y2 - p1_y)))
+    return 1;
+  return trn.complete_simple_poly();
+}
+
+
+int Item_func_buffer::Transporter::add_last_edge_buffer()
+{
+  Gcalc_operation_transporter trn(m_fn, m_heap);
+  double e1_x, e1_y, p1_x, p1_y;
+
+  ++m_nshapes;
+  if (trn.start_simple_poly())
+    return 1;
+
+  calculate_perpendicular(x1, y1, x2, y2, m_d, &e1_x, &e1_y, &p1_x, &p1_y);
+
+  if (trn.add_point(x1 + p1_x, y1 + p1_y) ||
+      trn.add_point(x1 - p1_x, y1 - p1_y) ||
+      trn.add_point(x2 - p1_x, y2 - p1_y) ||
+      fill_half_circle(&trn, x2, y2, -p1_x, -p1_y) ||
+      trn.add_point(x2 + p1_x, y2 + p1_y))
+    return 1;
+  return trn.complete_simple_poly();
+}
+
+
+int Item_func_buffer::Transporter::add_point_buffer(double x, double y)
+{
+  Gcalc_operation_transporter trn(m_fn, m_heap);
+
+  m_nshapes++;
+  if (trn.start_simple_poly())
+    return 1;
+  if (trn.add_point(x - m_d, y) ||
+      fill_half_circle(&trn, x, y, -m_d, 0.0) ||
+      trn.add_point(x + m_d, y) ||
+      fill_half_circle(&trn, x, y, m_d, 0.0))
+    return 1;
+  return trn.complete_simple_poly();
+}
+
+
+int Item_func_buffer::Transporter::start_line()
+{
+  m_npoints= 0;
+  int_start_line();
+  return 0;
+}
+
+
+int Item_func_buffer::Transporter::start_poly()
+{
+  ++m_nshapes;
+  return Gcalc_operation_transporter::start_poly();
+}
+
+
+int Item_func_buffer::Transporter::start_ring()
+{
+  m_npoints= 0;
+  return Gcalc_operation_transporter::start_ring();
+}
+
+
+int Item_func_buffer::Transporter::add_point(double x, double y)
+{
+  if (m_npoints && x == x2 && y == y2)
+    return 0;
+
+  ++m_npoints;
+
+  if (m_npoints == 1)
+  {
+    x00= x;
+    y00= y;
+  }
+  else if (m_npoints == 2)
+  {
+    x01= x;
+    y01= y;
+  }
+  else if (add_edge_buffer(x, y, (m_npoints == 3) && line_started(), false))
+    return 1;
+
+  x1= x2;
+  y1= y2;
+  x2= x;
+  y2= y;
+
+  return line_started() ? 0 : Gcalc_operation_transporter::add_point(x, y);
+}
+
+
+int Item_func_buffer::Transporter::complete()
+{
+  if (m_npoints)
+  {
+    if (m_npoints == 1)
+    {
+      if (add_point_buffer(x2, y2))
+        return 1;
+    }
+    else if (m_npoints == 2)
+    {
+      if (add_edge_buffer(x1, y1, true, true))
+        return 1;
+    }
+    else if (line_started())
+    {
+      if (add_last_edge_buffer())
+        return 1;
+    }
+    else
+    {
+      if (x2 != x00 && y2 != y00)
+      {
+        if (add_edge_buffer(x00, y00, false, false))
+          return 1;
+        x1= x2;
+        y1= y2;
+        x2= x00;
+        y2= y00;
+      }
+      if (add_edge_buffer(x01, y01, false, false))
+        return 1;
+    }
+  }
+
+  return 0;
+}
+
+
+int Item_func_buffer::Transporter::complete_line()
+{
+  if (complete())
+    return 1;
+  int_complete_line();
+  return 0;
+}
+
+
+int Item_func_buffer::Transporter::complete_ring()
+{
+  return complete() ||
+         Gcalc_operation_transporter::complete_ring();
+}
+
+
+String *Item_func_buffer::val_str(String *str_value)
+{
+  DBUG_ENTER("Item_func_buffer::val_str");
+  DBUG_ASSERT(fixed == 1);
+  String *obj= args[0]->val_str(&tmp_value);
+  double dist= args[1]->val_real();
+  Geometry_buffer buffer;
+  Geometry *g;
+  uint32 union_pos;
+  uint32 srid= 0;
+  String *str_result= NULL;
+  Transporter trn(&func, &collector, dist);
+
+  null_value= 1;
+  if (args[0]->null_value || args[1]->null_value ||
+      !(g= Geometry::construct(&buffer, obj->ptr(), obj->length())))
+    goto mem_error;
+
+  if (func.reserve_op_buffer(2))
+    goto mem_error;
+  /* will specify operands later */
+  union_pos= func.get_next_operation_pos();
+  func.add_operation((dist > 0.0) ? Gcalc_function::op_union :
+                                    Gcalc_function::op_difference, 0);
+
+  if (g->store_shapes(&trn))
+    goto mem_error;
+
+  func.add_operands_to_op(union_pos, trn.m_nshapes);
+  collector.prepare_operation();
+  if (func.alloc_states())
+    goto mem_error;
+  operation.init(&func);
+
+  if (operation.count_all(&collector) ||
+      operation.get_result(&res_receiver))
+    goto mem_error;
+
+
+  str_value->set_charset(&my_charset_bin);
+  if (str_value->reserve(SRID_SIZE, 512))
+    goto mem_error;
+  str_value->length(0);
+  str_value->q_append(srid);
+
+  if (!Geometry::create_from_opresult(&buffer, str_value, res_receiver))
+    goto mem_error;
+
+  null_value= 0;
+  str_result= str_value;
+mem_error:
+  collector.reset();
+  func.reset();
+  scan_it.reset();
+  res_receiver.reset();
+  DBUG_RETURN(str_result);
+}
+
+
 longlong Item_func_isempty::val_int()
 {
   DBUG_ASSERT(fixed == 1);
@@ -566,20 +1460,50 @@ longlong Item_func_isempty::val_int()
 }
 
 
-/**
-  @todo
-    Ramil or Holyfoot, add real IsSimple calculation
-*/
 longlong Item_func_issimple::val_int()
 {
-  DBUG_ASSERT(fixed == 1);
-  String tmp;
   String *swkb= args[0]->val_str(&tmp);
   Geometry_buffer buffer;
+  Gcalc_operation_transporter trn(&func, &collector);
+  Geometry *g;
+  int result= 1;
+
+  DBUG_ENTER("Item_func_issimple::val_int");
+  DBUG_ASSERT(fixed == 1);
   
-  null_value= args[0]->null_value ||
-              !(Geometry::construct(&buffer, swkb->ptr(), swkb->length()));
-  /* TODO: Ramil or Holyfoot, add real IsSimple calculation */
+  if ((null_value= args[0]->null_value) ||
+      !(g= Geometry::construct(&buffer, swkb->ptr(), swkb->length())))
+    DBUG_RETURN(0);
+
+
+  if (g->get_class_info()->m_type_id == Geometry::wkb_point)
+    DBUG_RETURN(1);
+
+  if (g->store_shapes(&trn))
+    goto mem_error;
+
+  collector.prepare_operation();
+  scan_it.init(&collector);
+
+  while (scan_it.more_points())
+  {
+    if (scan_it.step())
+      goto mem_error;
+
+    if (scan_it.get_event() == scev_intersection)
+    {
+      result= 0;
+      break;
+    }
+  }
+
+  collector.reset();
+  func.reset();
+  scan_it.reset();
+  DBUG_RETURN(result);
+mem_error:
+  null_value= 1;
+  DBUG_RETURN(0);
   return 0;
 }
 
@@ -752,4 +1676,160 @@ longlong Item_func_srid::val_int()
   return (longlong) (uint4korr(swkb->ptr()));
 }
 
+
+double Item_func_distance::val_real()
+{
+  bool above_cur_point, cur_point_edge;
+  const Gcalc_scan_iterator::point *evpos;
+  const Gcalc_heap::Info *cur_point, *dist_point;
+  Gcalc_scan_events ev;
+  double t, distance, cur_distance;
+  double x1, x2, y1, y2;
+  double ex, ey, vx, vy, e_sqrlen;
+  uint obj2_si;
+  Gcalc_operation_transporter trn(&func, &collector);
+
+  DBUG_ENTER("Item_func_distance::val_real");
+  DBUG_ASSERT(fixed == 1);
+  String *res1= args[0]->val_str(&tmp_value1);
+  String *res2= args[1]->val_str(&tmp_value2);
+  Geometry_buffer buffer1, buffer2;
+  Geometry *g1, *g2;
+
+
+  if ((null_value= (args[0]->null_value || args[1]->null_value ||
+          !(g1= Geometry::construct(&buffer1, res1->ptr(), res1->length())) ||
+          !(g2= Geometry::construct(&buffer2, res2->ptr(), res2->length())))))
+    goto mem_error;
+
+  if ((g1->get_class_info()->m_type_id == Geometry::wkb_point) &&
+      (g2->get_class_info()->m_type_id == Geometry::wkb_point))
+  {
+    if (((Gis_point *) g1)->get_xy(&x1, &y1) ||
+        ((Gis_point *) g2)->get_xy(&x2, &y2))
+      goto mem_error;
+    ex= x2 - x1;
+    ey= y2 - y1;
+    DBUG_RETURN(sqrt(ex * ex + ey * ey));
+  }
+
+  if (func.reserve_op_buffer(1))
+    goto mem_error;
+  func.add_operation(Gcalc_function::op_intersection, 2);
+
+  if (g1->store_shapes(&trn))
+    goto mem_error;
+  obj2_si= func.get_nshapes();
+  if (g2->store_shapes(&trn) || func.alloc_states())
+    goto mem_error;
+
+  collector.prepare_operation();
+  scan_it.init(&collector);
+
+  above_cur_point= false;
+  distance= DBL_MAX;
+  while (scan_it.more_points())
+  {
+    if (scan_it.step())
+      goto mem_error;
+    evpos= scan_it.get_event_position();
+    ev= scan_it.get_event();
+    cur_point= evpos->pi;
+
+    /*
+       handling intersection we only need to check if it's the intersecion
+       of objects 1 and 2. In this case distance is 0
+    */
+    if (ev == scev_intersection)
+    {
+      if ((evpos->get_next()->pi->shape >= obj2_si) !=
+            (cur_point->shape >= obj2_si))
+      {
+        distance= 0;
+        goto exit;
+      }
+      continue;
+    }
+
+    /*
+       if we get 'scev_point | scev_end | scev_two_ends' we don't need
+       to check for intersection of objects.
+       Though we need to calculate distances.
+    */
+    if (ev & (scev_point | scev_end | scev_two_ends))
+      goto count_distance;
+
+    /*
+       having these events we need to check for possible intersection
+       of objects
+       scev_thread | scev_two_threads | scev_single_point
+    */
+    DBUG_ASSERT(ev & (scev_thread | scev_two_threads | scev_single_point));
+
+    func.clear_state();
+    for (Gcalc_point_iterator pit(&scan_it); pit.point() != evpos; ++pit)
+    {
+      gcalc_shape_info si= pit.point()->get_shape();
+      if ((func.get_shape_kind(si) == Gcalc_function::shape_polygon))
+        func.invert_state(si);
+    }
+    func.invert_state(evpos->get_shape());
+    if (func.count())
+    {
+      /* Point of one object is inside the other - intersection found */
+      distance= 0;
+      goto exit;
+    }
+
+
+count_distance:
+    if (cur_point->shape >= obj2_si)
+      continue;
+    cur_point_edge= !cur_point->is_bottom();
+
+    for (dist_point= collector.get_first(); dist_point; dist_point= dist_point->get_next())
+    {
+      /* We only check vertices of object 2 */
+      if (dist_point->shape < obj2_si)
+        continue;
+
+      /* if we have an edge to check */
+      if (dist_point->left)
+      {
+        t= count_edge_t(dist_point, dist_point->left, cur_point,
+                        ex, ey, vx, vy, e_sqrlen);
+        if ((t>0.0) && (t<1.0))
+        {
+          cur_distance= distance_to_line(ex, ey, vx, vy, e_sqrlen);
+          if (distance > cur_distance)
+            distance= cur_distance;
+        }
+      }
+      if (cur_point_edge)
+      {
+        t= count_edge_t(cur_point, cur_point->left, dist_point,
+                        ex, ey, vx, vy, e_sqrlen);
+        if ((t>0.0) && (t<1.0))
+        {
+          cur_distance= distance_to_line(ex, ey, vx, vy, e_sqrlen);
+          if (distance > cur_distance)
+            distance= cur_distance;
+        }
+      }
+      cur_distance= distance_points(cur_point, dist_point);
+      if (distance > cur_distance)
+        distance= cur_distance;
+    }
+  }
+exit:
+  collector.reset();
+  func.reset();
+  scan_it.reset();
+  DBUG_RETURN(distance);
+mem_error:
+  null_value= 1;
+  DBUG_RETURN(0);
+}
+
+
 #endif /*HAVE_SPATIAL*/
diff --git a/sql/item_geofunc.h b/sql/item_geofunc.h
index 0bcd933e52b..1838a3783e8 100644
--- a/sql/item_geofunc.h
+++ b/sql/item_geofunc.h
@@ -1,4 +1,7 @@
-/* Copyright (c) 2003, 2011, Oracle and/or its affiliates.
+#ifndef ITEM_GEOFUNC_INCLUDED
+#define ITEM_GEOFUNC_INCLUDED
+
+/* Copyright (c) 2000, 2010, Oracle and/or its affiliates. All rights reserved.
 
    This program is free software; you can redistribute it and/or modify
    it under the terms of the GNU General Public License as published by
@@ -10,8 +13,8 @@
    GNU General Public License for more details.
 
    You should have received a copy of the GNU General Public License
-   along with this program; if not, write to the Free Software
-   Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA */
+   along with this program; if not, write to the Free Software Foundation,
+   51 Franklin Street, Suite 500, Boston, MA 02110-1335 USA */
 
 
 /* This file defines all spatial functions */
@@ -22,6 +25,8 @@
 #pragma interface			/* gcc class implementation */
 #endif
 
+#include "gcalc_slicescan.h"
+
 class Item_geometry_func: public Item_str_func
 {
 public:
@@ -41,7 +46,7 @@ class Item_func_geometry_from_text: public Item_geometry_func
 public:
   Item_func_geometry_from_text(Item *a) :Item_geometry_func(a) {}
   Item_func_geometry_from_text(Item *a, Item *srid) :Item_geometry_func(a, srid) {}
-  const char *func_name() const { return "geometryfromtext"; }
+  const char *func_name() const { return "st_geometryfromtext"; }
   String *val_str(String *);
 };
 
@@ -50,7 +55,7 @@ class Item_func_geometry_from_wkb: public Item_geometry_func
 public:
   Item_func_geometry_from_wkb(Item *a): Item_geometry_func(a) {}
   Item_func_geometry_from_wkb(Item *a, Item *srid): Item_geometry_func(a, srid) {}
-  const char *func_name() const { return "geometryfromwkb"; }
+  const char *func_name() const { return "st_geometryfromwkb"; }
   String *val_str(String *);
 };
 
@@ -58,7 +63,7 @@ class Item_func_as_wkt: public Item_str_func
 {
 public:
   Item_func_as_wkt(Item *a): Item_str_func(a) {}
-  const char *func_name() const { return "astext"; }
+  const char *func_name() const { return "st_astext"; }
   String *val_str(String *);
   void fix_length_and_dec();
 };
@@ -67,7 +72,7 @@ class Item_func_as_wkb: public Item_geometry_func
 {
 public:
   Item_func_as_wkb(Item *a): Item_geometry_func(a) {}
-  const char *func_name() const { return "aswkb"; }
+  const char *func_name() const { return "st_aswkb"; }
   String *val_str(String *);
   enum_field_types field_type() const  { return MYSQL_TYPE_BLOB; }
 };
@@ -77,11 +82,12 @@ class Item_func_geometry_type: public Item_str_func
 public:
   Item_func_geometry_type(Item *a): Item_str_func(a) {}
   String *val_str(String *);
-  const char *func_name() const { return "geometrytype"; }
+  const char *func_name() const { return "st_geometrytype"; }
   void fix_length_and_dec() 
   {
-     max_length=20; // "GeometryCollection" is the most long
-     maybe_null= 1;
+    // "GeometryCollection" is the longest
+    max_length= 20;
+    maybe_null= 1;
   };
 };
 
@@ -89,7 +95,7 @@ class Item_func_centroid: public Item_geometry_func
 {
 public:
   Item_func_centroid(Item *a): Item_geometry_func(a) {}
-  const char *func_name() const { return "centroid"; }
+  const char *func_name() const { return "st_centroid"; }
   String *val_str(String *);
   Field::geometry_type get_geometry_type() const;
 };
@@ -98,7 +104,7 @@ class Item_func_envelope: public Item_geometry_func
 {
 public:
   Item_func_envelope(Item *a): Item_geometry_func(a) {}
-  const char *func_name() const { return "envelope"; }
+  const char *func_name() const { return "st_envelope"; }
   String *val_str(String *);
   Field::geometry_type get_geometry_type() const;
 };
@@ -108,7 +114,7 @@ class Item_func_point: public Item_geometry_func
 public:
   Item_func_point(Item *a, Item *b): Item_geometry_func(a, b) {}
   Item_func_point(Item *a, Item *b, Item *srid): Item_geometry_func(a, b, srid) {}
-  const char *func_name() const { return "point"; }
+  const char *func_name() const { return "st_point"; }
   String *val_str(String *);
   Field::geometry_type get_geometry_type() const;
 };
@@ -124,11 +130,11 @@ public:
     switch (decomp_func)
     {
       case SP_STARTPOINT:
-        return "startpoint";
+        return "st_startpoint";
       case SP_ENDPOINT:
-        return "endpoint";
+        return "st_endpoint";
       case SP_EXTERIORRING:
-        return "exteriorring";
+        return "st_exteriorring";
       default:
 	DBUG_ASSERT(0);  // Should never happened
         return "spatial_decomp_unknown"; 
@@ -148,11 +154,11 @@ public:
     switch (decomp_func_n)
     {
       case SP_POINTN:
-        return "pointn";
+        return "st_pointn";
       case SP_GEOMETRYN:
-        return "geometryn";
+        return "st_geometryn";
       case SP_INTERIORRINGN:
-        return "interiorringn";
+        return "st_interiorringn";
       default:
 	DBUG_ASSERT(0);  // Should never happened
         return "spatial_decomp_n_unknown"; 
@@ -177,7 +183,6 @@ public:
   String *val_str(String *);
   void fix_length_and_dec()
   {
-    Item_geometry_func::fix_length_and_dec();
     for (unsigned int i= 0; i < arg_count; ++i)
     {
       if (args[i]->fixed && args[i]->field_type() != MYSQL_TYPE_GEOMETRY)
@@ -191,57 +196,53 @@ public:
     }
   }
  
-  const char *func_name() const { return "multipoint"; }
+  const char *func_name() const { return "st_multipoint"; }
 };
 
+
 /*
   Spatial relations
 */
 
-class Item_func_spatial_rel: public Item_bool_func2
+class Item_func_spatial_mbr_rel: public Item_bool_func2
 {
   enum Functype spatial_rel;
 public:
-  Item_func_spatial_rel(Item *a,Item *b, enum Functype sp_rel) :
+  Item_func_spatial_mbr_rel(Item *a,Item *b, enum Functype sp_rel) :
     Item_bool_func2(a,b) { spatial_rel = sp_rel; }
   longlong val_int();
   enum Functype functype() const 
   { 
-    switch (spatial_rel) {
-    case SP_CONTAINS_FUNC:
-      return SP_WITHIN_FUNC;
-    case SP_WITHIN_FUNC:
-      return SP_CONTAINS_FUNC;
-    default:
-      return spatial_rel;
-    }
+    return spatial_rel;
   }
   enum Functype rev_functype() const { return spatial_rel; }
-  const char *func_name() const 
-  { 
-    switch (spatial_rel) {
-    case SP_CONTAINS_FUNC:
-      return "contains";
-    case SP_WITHIN_FUNC:
-      return "within";
-    case SP_EQUALS_FUNC:
-      return "equals";
-    case SP_DISJOINT_FUNC:
-      return "disjoint";
-    case SP_INTERSECTS_FUNC:
-      return "intersects";
-    case SP_TOUCHES_FUNC:
-      return "touches";
-    case SP_CROSSES_FUNC:
-      return "crosses";
-    case SP_OVERLAPS_FUNC:
-      return "overlaps";
-    default:
-      DBUG_ASSERT(0);  // Should never happened
-      return "sp_unknown"; 
-    }
+  const char *func_name() const;
+  virtual inline void print(String *str, enum_query_type query_type)
+  {
+    Item_func::print(str, query_type);
   }
+  void fix_length_and_dec() { maybe_null= 1; }
+  bool is_null() { (void) val_int(); return null_value; }
+};
+
 
+class Item_func_spatial_rel: public Item_bool_func2
+{
+  enum Functype spatial_rel;
+  Gcalc_heap collector;
+  Gcalc_scan_iterator scan_it;
+  Gcalc_function func;
+  String tmp_value1,tmp_value2;
+public:
+  Item_func_spatial_rel(Item *a,Item *b, enum Functype sp_rel);
+  virtual ~Item_func_spatial_rel();
+  longlong val_int();
+  enum Functype functype() const 
+  { 
+    return spatial_rel;
+  }
+  enum Functype rev_functype() const { return spatial_rel; }
+  const char *func_name() const;
   virtual inline void print(String *str, enum_query_type query_type)
   {
     Item_func::print(str, query_type);
@@ -249,25 +250,104 @@ public:
 
   void fix_length_and_dec() { maybe_null= 1; }
   bool is_null() { (void) val_int(); return null_value; }
+protected:
+  int func_touches();
+  int func_equals();
 };
 
+
+/*
+  Spatial operations
+*/
+
+class Item_func_spatial_operation: public Item_geometry_func
+{
+public:
+  Gcalc_function::op_type spatial_op;
+  Gcalc_heap collector;
+  Gcalc_function func;
+  Gcalc_scan_iterator scan_it;
+
+  Gcalc_result_receiver res_receiver;
+  Gcalc_operation_reducer operation;
+  String tmp_value1,tmp_value2;
+public:
+  Item_func_spatial_operation(Item *a,Item *b, Gcalc_function::op_type sp_op) :
+    Item_geometry_func(a, b), spatial_op(sp_op)
+  {}
+  virtual ~Item_func_spatial_operation();
+  String *val_str(String *);
+  const char *func_name() const;
+  virtual inline void print(String *str, enum_query_type query_type)
+  {
+    Item_func::print(str, query_type);
+  }
+};
+
+
+class Item_func_buffer: public Item_geometry_func
+{
+protected:
+  class Transporter : public Gcalc_operation_transporter
+  {
+    int m_npoints;
+    double m_d;
+    double x1,y1,x2,y2;
+    double x00,y00,x01,y01;
+    int add_edge_buffer(double x3, double y3, bool round_p1, bool round_p2);
+    int add_last_edge_buffer();
+    int add_point_buffer(double x, double y);
+    int complete();
+  public:
+    int m_nshapes;
+    Transporter(Gcalc_function *fn, Gcalc_heap *heap, double d) :
+      Gcalc_operation_transporter(fn, heap), m_npoints(0), m_d(d), m_nshapes(0)
+    {}
+    int single_point(double x, double y);
+    int start_line();
+    int complete_line();
+    int start_poly();
+    int start_ring();
+    int complete_ring();
+    int add_point(double x, double y);
+  };
+  Gcalc_heap collector;
+  Gcalc_function func;
+  Gcalc_scan_iterator scan_it;
+
+  Gcalc_result_receiver res_receiver;
+  Gcalc_operation_reducer operation;
+  String tmp_value;
+
+public:
+  Item_func_buffer(Item *obj, Item *distance):
+    Item_geometry_func(obj, distance) {}
+  const char *func_name() const { return "st_buffer"; }
+  String *val_str(String *);
+};
+
+
 class Item_func_isempty: public Item_bool_func
 {
 public:
   Item_func_isempty(Item *a): Item_bool_func(a) {}
   longlong val_int();
   optimize_type select_optimize() const { return OPTIMIZE_NONE; }
-  const char *func_name() const { return "isempty"; }
+  const char *func_name() const { return "st_isempty"; }
   void fix_length_and_dec() { maybe_null= 1; }
 };
 
 class Item_func_issimple: public Item_bool_func
 {
+  Gcalc_heap collector;
+  Gcalc_function func;
+  Gcalc_scan_iterator scan_it;
+  String tmp;
 public:
   Item_func_issimple(Item *a): Item_bool_func(a) {}
   longlong val_int();
   optimize_type select_optimize() const { return OPTIMIZE_NONE; }
-  const char *func_name() const { return "issimple"; }
+  const char *func_name() const { return "st_issimple"; }
   void fix_length_and_dec() { maybe_null= 1; }
 };
 
@@ -277,7 +357,7 @@ public:
   Item_func_isclosed(Item *a): Item_bool_func(a) {}
   longlong val_int();
   optimize_type select_optimize() const { return OPTIMIZE_NONE; }
-  const char *func_name() const { return "isclosed"; }
+  const char *func_name() const { return "st_isclosed"; }
   void fix_length_and_dec() { maybe_null= 1; }
 };
 
@@ -287,7 +367,7 @@ class Item_func_dimension: public Item_int_func
 public:
   Item_func_dimension(Item *a): Item_int_func(a) {}
   longlong val_int();
-  const char *func_name() const { return "dimension"; }
+  const char *func_name() const { return "st_dimension"; }
   void fix_length_and_dec() { max_length= 10; maybe_null= 1; }
 };
 
@@ -297,7 +377,7 @@ class Item_func_x: public Item_real_func
 public:
   Item_func_x(Item *a): Item_real_func(a) {}
   double val_real();
-  const char *func_name() const { return "x"; }
+  const char *func_name() const { return "st_x"; }
   void fix_length_and_dec() 
   { 
     Item_real_func::fix_length_and_dec();
@@ -312,7 +392,7 @@ class Item_func_y: public Item_real_func
 public:
   Item_func_y(Item *a): Item_real_func(a) {}
   double val_real();
-  const char *func_name() const { return "y"; }
+  const char *func_name() const { return "st_y"; }
   void fix_length_and_dec() 
   { 
     Item_real_func::fix_length_and_dec();
@@ -327,7 +407,7 @@ class Item_func_numgeometries: public Item_int_func
 public:
   Item_func_numgeometries(Item *a): Item_int_func(a) {}
   longlong val_int();
-  const char *func_name() const { return "numgeometries"; }
+  const char *func_name() const { return "st_numgeometries"; }
   void fix_length_and_dec() { max_length= 10; maybe_null= 1; }
 };
 
@@ -338,7 +418,7 @@ class Item_func_numinteriorring: public Item_int_func
 public:
   Item_func_numinteriorring(Item *a): Item_int_func(a) {}
   longlong val_int();
-  const char *func_name() const { return "numinteriorrings"; }
+  const char *func_name() const { return "st_numinteriorrings"; }
   void fix_length_and_dec() { max_length= 10; maybe_null= 1; }
 };
 
@@ -349,7 +429,7 @@ class Item_func_numpoints: public Item_int_func
 public:
   Item_func_numpoints(Item *a): Item_int_func(a) {}
   longlong val_int();
-  const char *func_name() const { return "numpoints"; }
+  const char *func_name() const { return "st_numpoints"; }
   void fix_length_and_dec() { max_length= 10; maybe_null= 1; }
 };
 
@@ -360,7 +440,7 @@ class Item_func_area: public Item_real_func
 public:
   Item_func_area(Item *a): Item_real_func(a) {}
   double val_real();
-  const char *func_name() const { return "area"; }
+  const char *func_name() const { return "st_area"; }
   void fix_length_and_dec() 
   { 
     Item_real_func::fix_length_and_dec();
@@ -375,7 +455,7 @@ class Item_func_glength: public Item_real_func
 public:
   Item_func_glength(Item *a): Item_real_func(a) {}
   double val_real();
-  const char *func_name() const { return "glength"; }
+  const char *func_name() const { return "st_length"; }
   void fix_length_and_dec() 
   { 
     Item_real_func::fix_length_and_dec();
@@ -394,11 +474,26 @@ public:
   void fix_length_and_dec() { max_length= 10; maybe_null= 1; }
 };
 
+
+class Item_func_distance: public Item_real_func
+{
+  String tmp_value1;
+  String tmp_value2;
+  Gcalc_heap collector;
+  Gcalc_function func;
+  Gcalc_scan_iterator scan_it;
+public:
+  Item_func_distance(Item *a, Item *b): Item_real_func(a, b) {}
+  double val_real();
+  const char *func_name() const { return "st_distance"; }
+};
+
 #define GEOM_NEW(thd, obj_constructor) new (thd->mem_root) obj_constructor
 
 #else /*HAVE_SPATIAL*/
 
 #define GEOM_NEW(thd, obj_constructor) NULL
 
-#endif
+#endif /*HAVE_SPATIAL*/
+#endif /*ITEM_GEOFUNC_INCLUDED*/
 
diff --git a/sql/plistsort.c b/sql/plistsort.c
new file mode 100644
index 00000000000..71d287e7b45
--- /dev/null
+++ b/sql/plistsort.c
@@ -0,0 +1,166 @@
+/* Copyright (c) 2000, 2010 Oracle and/or its affiliates. All rights reserved.
+
+   This program is free software; you can redistribute it and/or modify
+   it under the terms of the GNU General Public License as published by
+   the Free Software Foundation; version 2 of the License.
+
+   This program is distributed in the hope that it will be useful,
+   but WITHOUT ANY WARRANTY; without even the implied warranty of
+   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+   GNU General Public License for more details.
+
+   You should have received a copy of the GNU General Public License
+   along with this program; if not, write to the Free Software
+   Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA */
+
+
+/*
+things to define before including the file:
+
+#define LS_LIST_ITEM ListItem
+#define LS_COMPARE_FUNC_DECL compare_func var_name,
+#define LS_COMPARE_FUNC_CALL(list_el1, list_el2) (*var_name)(list_el1, list_el2)
+#define LS_NEXT(A) (A)->next
+#define LS_SET_NEXT(A,val) (A)->next= val
+#define LS_P_NEXT(A) &(A)->next
+#define LS_NAME plistsort
+#define LS_SCOPE static
+#define LS_STRUCT_NAME ls_struct_name
+*/
+
+typedef struct LS_STRUCT_NAME
+{
+  LS_LIST_ITEM *list1;
+  int list_len;
+  int return_point;
+} LS_STRUCT_NAME;
+
+LS_SCOPE LS_LIST_ITEM* LS_NAME(LS_COMPARE_FUNC_DECL LS_LIST_ITEM *list, int list_len)
+{
+  LS_LIST_ITEM *list_end;
+  LS_LIST_ITEM *sorted_list;
+
+  struct LS_STRUCT_NAME stack[63], *sp= stack;
+
+  if (list_len < 2)
+    return list;
+
+  sp->list_len= list_len;
+  sp->return_point= 2;
+
+recursion_point:
+
+  if (sp->list_len < 4)
+  {
+    LS_LIST_ITEM *e1, *e2;
+    sorted_list= list;
+    e1= LS_NEXT(sorted_list);
+    list_end= LS_NEXT(e1);
+    if (LS_COMPARE_FUNC_CALL(sorted_list, e1))
+    {
+      sorted_list= e1;
+      e1= list;
+    }
+    if (sp->list_len == 2)
+    {
+      LS_SET_NEXT(sorted_list, e1);
+      LS_SET_NEXT(e1, NULL);
+      goto exit_point;
+    }
+    e2= list_end;
+    list_end= LS_NEXT(e2);
+    if (LS_COMPARE_FUNC_CALL(e1, e2))
+    {
+      {
+        LS_LIST_ITEM *tmp_e= e1;
+        e1= e2;
+        e2= tmp_e;
+      }
+      if (LS_COMPARE_FUNC_CALL(sorted_list, e1))
+      {
+        LS_LIST_ITEM *tmp_e= sorted_list;
+        sorted_list= e1;
+        e1= tmp_e;
+      }
+    }
+
+    LS_SET_NEXT(sorted_list, e1);
+    LS_SET_NEXT(e1, e2);
+    LS_SET_NEXT(e2, NULL);
+    goto exit_point;
+  }
+
+  {
+    register struct LS_STRUCT_NAME *sp0= sp++;
+    sp->list_len= sp0->list_len >> 1;
+    sp0->list_len-= sp->list_len;
+    sp->return_point= 0;
+  }
+  goto recursion_point;
+return_point0:
+  sp->list1= sorted_list;
+  {
+    register struct LS_STRUCT_NAME *sp0= sp++;
+    list= list_end;
+    sp->list_len= sp0->list_len;
+    sp->return_point= 1;
+  }
+  goto recursion_point;
+return_point1:
+  {
+    register LS_LIST_ITEM **hook= &sorted_list;
+    register LS_LIST_ITEM *list1= sp->list1;
+    register LS_LIST_ITEM *list2= sorted_list;
+
+    if (LS_COMPARE_FUNC_CALL(list1, list2))
+    {
+      LS_LIST_ITEM *tmp_e= list2;
+      list2= list1;
+      list1= tmp_e;
+    }
+    for (;;)
+    {
+      *hook= list1;
+      do
+      {
+        if (!(list1= *(hook= LS_P_NEXT(list1))))
+        {
+          *hook= list2;
+          goto exit_point;
+        }
+      } while (LS_COMPARE_FUNC_CALL(list2, list1));
+
+      *hook= list2;
+      do
+      {
+        if (!(list2= *(hook= LS_P_NEXT(list2))))
+        {
+          *hook= list1;
+          goto exit_point;
+        }
+      } while (LS_COMPARE_FUNC_CALL(list1, list2));
+    }
+  }
+
+exit_point:
+  switch ((sp--)->return_point)
+  {
+    case 0: goto return_point0;
+    case 1: goto return_point1;
+    default:;
+  }
+
+  return sorted_list;
+}
+
+
+#undef LS_LIST_ITEM
+#undef LS_NEXT
+#undef LS_SET_NEXT
+#undef LS_P_NEXT
+#undef LS_NAME
+#undef LS_STRUCT_NAME
+#undef LS_SCOPE
+#undef LS_COMPARE_FUNC_DECL
+#undef LS_COMPARE_FUNC_CALL
+
diff --git a/sql/spatial.cc b/sql/spatial.cc
index 8b869a5b1ca..48f6f466a30 100644
--- a/sql/spatial.cc
+++ b/sql/spatial.cc
@@ -1,4 +1,4 @@
-/* Copyright (C) 2004 MySQL AB
+/* Copyright (c) 2004, 2010, Oracle and/or its affiliates. All rights reserved.
 
    This program is free software; you can redistribute it and/or modify
    it under the terms of the GNU General Public License as published by
@@ -10,8 +10,8 @@
    GNU General Public License for more details.
 
    You should have received a copy of the GNU General Public License
-   along with this program; if not, write to the Free Software
-   Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA */
+   along with this program; if not, write to the Free Software Foundation,
+   51 Franklin Street, Suite 500, Boston, MA 02110-1335 USA */
 
 #include "mysql_priv.h"
 
@@ -38,7 +38,7 @@
   17
 */
 
-#define MAX_DIGITS_IN_DOUBLE 22
+#define MAX_DIGITS_IN_DOUBLE 30
 
 /***************************** Gis_class_info *******************************/
 
@@ -264,12 +264,28 @@ Geometry *Geometry::create_from_wkb(Geometry_buffer *buffer,
 }
 
 
+int Geometry::create_from_opresult(Geometry_buffer *g_buf,
+                                   String *res, Gcalc_result_receiver &rr)
+{
+  uint32 geom_type= rr.get_result_typeid();
+  Geometry *obj= create_by_typeid(g_buf, geom_type);
+
+  if (!obj || res->reserve(WKB_HEADER_SIZE, 512))
+    return 1;
+
+  res->q_append((char) wkb_ndr);
+  res->q_append(geom_type);
+  return obj->init_from_opresult(res, rr.result(), rr.get_nshapes());
+}
+
+
 bool Geometry::envelope(String *result) const
 {
   MBR mbr;
   const char *end;
 
-  if (get_mbr(&mbr, &end) || result->reserve(1+4*3+SIZEOF_STORED_DOUBLE*10))
+  if (get_mbr(&mbr, &end) ||
+      result->reserve(1 + 4 * 3 + SIZEOF_STORED_DOUBLE * 10))
     return 1;
 
   result->q_append((char) wkb_ndr);
@@ -306,13 +322,13 @@ bool Geometry::envelope(String *result) const
 
 bool Geometry::create_point(String *result, const char *data) const
 {
-  if (no_data(data, SIZEOF_STORED_DOUBLE * 2) ||
-      result->reserve(1 + 4 + SIZEOF_STORED_DOUBLE * 2))
+  if (no_data(data, POINT_DATA_SIZE) ||
+      result->reserve(1 + 4 + POINT_DATA_SIZE))
     return 1;
   result->q_append((char) wkb_ndr);
   result->q_append((uint32) wkb_point);
   /* Copy two double in same format */
-  result->q_append(data, SIZEOF_STORED_DOUBLE*2);
+  result->q_append(data, POINT_DATA_SIZE);
   return 0;
 }
 
@@ -332,7 +348,7 @@ bool Geometry::create_point(String *result, const char *data) const
 
 bool Geometry::create_point(String *result, double x, double y) const
 {
-  if (result->reserve(1 + 4 + SIZEOF_STORED_DOUBLE * 2))
+  if (result->reserve(1 + 4 + POINT_DATA_SIZE))
     return 1;
 
   result->q_append((char) wkb_ndr);
@@ -364,7 +380,7 @@ const char *Geometry::append_points(String *txt, uint32 n_points,
     double x,y;
     data+= offset;
     get_point(&x, &y, data);
-    data+= SIZEOF_STORED_DOUBLE * 2;
+    data+= POINT_DATA_SIZE;
     txt->qs_append(x);
     txt->qs_append(' ');
     txt->qs_append(y);
@@ -399,7 +415,7 @@ const char *Geometry::get_mbr_for_points(MBR *mbr, const char *data,
   points= uint4korr(data);
   data+= 4;
 
-  if (no_data(data, (SIZEOF_STORED_DOUBLE * 2 + offset) * points))
+  if (no_data(data, (POINT_DATA_SIZE + offset) * points))
     return 0;
 
   /* Calculate MBR for points */
@@ -407,7 +423,7 @@ const char *Geometry::get_mbr_for_points(MBR *mbr, const char *data,
   {
     data+= offset;
     mbr->add_xy(data, data + SIZEOF_STORED_DOUBLE);
-    data+= SIZEOF_STORED_DOUBLE * 2;
+    data+= POINT_DATA_SIZE;
   }
   return data;
 }
@@ -425,7 +441,7 @@ bool Gis_point::init_from_wkt(Gis_read_stream *trs, String *wkb)
 {
   double x, y;
   if (trs->get_next_number(&x) || trs->get_next_number(&y) ||
-      wkb->reserve(SIZEOF_STORED_DOUBLE * 2))
+      wkb->reserve(POINT_DATA_SIZE))
     return 1;
   wkb->q_append(x);
   wkb->q_append(y);
@@ -472,6 +488,15 @@ bool Gis_point::get_mbr(MBR *mbr, const char **end) const
   return 0;
 }
 
+
+int Gis_point::store_shapes(Gcalc_shape_transporter *trn) const
+{
+  double x, y;
+
+  return get_xy(&x, &y) || trn->single_point(x, y);
+}
+
+
 const Geometry::Class_info *Gis_point::get_class_info() const
 {
   return &point_class;
@@ -523,12 +548,12 @@ uint Gis_line_string::init_from_wkb(const char *wkb, uint len,
   const char *wkb_end;
   Gis_point p;
 
-  if (len < 4)
+  if (len < 4 ||
+      (n_points= wkb_get_uint(wkb, bo))<1)
     return 0;
-  n_points= wkb_get_uint(wkb, bo);
   proper_length= 4 + n_points * POINT_DATA_SIZE;
 
-  if (!n_points || len < proper_length || res->reserve(proper_length))
+  if (len < proper_length || res->reserve(proper_length))
     return 0;
 
   res->q_append(n_points);
@@ -554,7 +579,7 @@ bool Gis_line_string::get_data_as_wkt(String *txt, const char **end) const
   data += 4;
 
   if (n_points < 1 ||
-      no_data(data, SIZEOF_STORED_DOUBLE * 2 * n_points) ||
+      no_data(data, POINT_DATA_SIZE * n_points) ||
       txt->reserve(((MAX_DIGITS_IN_DOUBLE + 1)*2 + 1) * n_points))
     return 1;
 
@@ -562,7 +587,7 @@ bool Gis_line_string::get_data_as_wkt(String *txt, const char **end) const
   {
     double x, y;
     get_point(&x, &y, data);
-    data+= SIZEOF_STORED_DOUBLE * 2;
+    data+= POINT_DATA_SIZE;
     txt->qs_append(x);
     txt->qs_append(' ');
     txt->qs_append(y);
@@ -591,17 +616,16 @@ int Gis_line_string::geom_length(double *len) const
     return 1;
   n_points= uint4korr(data);
   data+= 4;
-  if (n_points < 1 || no_data(data, SIZEOF_STORED_DOUBLE * 2 * n_points))
+  if (n_points < 1 || no_data(data, POINT_DATA_SIZE * n_points))
     return 1;
 
   get_point(&prev_x, &prev_y, data);
-  data+= SIZEOF_STORED_DOUBLE*2;
-
+  data+= POINT_DATA_SIZE;
   while (--n_points)
   {
     double x, y;
     get_point(&x, &y, data);
-    data+= SIZEOF_STORED_DOUBLE * 2;
+    data+= POINT_DATA_SIZE;
     *len+= sqrt(pow(prev_x-x,2)+pow(prev_y-y,2));
     prev_x= x;
     prev_y= y;
@@ -625,14 +649,14 @@ int Gis_line_string::is_closed(int *closed) const
     return 0;
   }
   data+= 4;
-  if (no_data(data, SIZEOF_STORED_DOUBLE * 2 * n_points))
+  if (no_data(data, POINT_DATA_SIZE * n_points))
     return 1;
 
   /* Get first point */
   get_point(&x1, &y1, data);
 
   /* get last point */
-  data+= SIZEOF_STORED_DOUBLE*2 + (n_points-2)*POINT_DATA_SIZE;
+  data+= POINT_DATA_SIZE + (n_points-2)*POINT_DATA_SIZE;
   get_point(&x2, &y2, data);
 
   *closed= (x1==x2) && (y1==y2);
@@ -676,6 +700,34 @@ int Gis_line_string::point_n(uint32 num, String *result) const
   return create_point(result, m_data + 4 + (num - 1) * POINT_DATA_SIZE);
 }
 
+
+int Gis_line_string::store_shapes(Gcalc_shape_transporter *trn) const
+{
+  uint32 n_points;
+  double x, y;
+  const char *data= m_data;
+
+  if (no_data(m_data, 4))
+    return 1;
+  n_points= uint4korr(data);
+  data+= 4;
+  if (n_points < 1 || no_data(data, POINT_DATA_SIZE * n_points))
+    return 1;
+
+  trn->start_line();
+
+  while (n_points--)
+  {
+    get_point(&x, &y, data);
+    data+= POINT_DATA_SIZE;
+    if (trn->add_point(x, y))
+      return 1;
+  }
+
+  return trn->complete_line();
+}
+
+
 const Geometry::Class_info *Gis_line_string::get_class_info() const
 {
   return &linestring_class;
@@ -737,6 +789,53 @@ bool Gis_polygon::init_from_wkt(Gis_read_stream *trs, String *wkb)
 }
 
 
+uint Gis_polygon::priv_init_from_opresult(String *bin,
+                                          const char *opres, uint32 n_shapes,
+                                          uint32 *poly_shapes)
+{
+  const char *opres_orig= opres;
+  uint32 position= bin->length();
+
+  *poly_shapes= 0;
+  if (bin->reserve(4, 512))
+    return 0;
+  bin->q_append(*poly_shapes);
+
+  while (n_shapes--)
+  {
+    uint32 n_points, proper_length;
+    const char *op_end, *p1_position;
+    Gis_point p;
+    Gcalc_function::shape_type st;
+
+    st= (Gcalc_function::shape_type) uint4korr(opres);
+    if (*poly_shapes && st != Gcalc_function::shape_hole)
+      break;
+    (*poly_shapes)++;
+    n_points= uint4korr(opres + 4) + 1; /* skip shape type id */
+    proper_length= 4 + n_points * POINT_DATA_SIZE;
+
+    if (bin->reserve(proper_length, 512))
+      return 0;
+
+    bin->q_append(n_points);
+    op_end= opres + 8 + (n_points-1) * 8 * 2;
+    p1_position= (opres+= 8);
+    for (; opres<op_end; opres+= POINT_DATA_SIZE)
+    {
+      if (!p.init_from_wkb(opres, POINT_DATA_SIZE, wkb_ndr, bin))
+        return 0;
+    }
+    if (!p.init_from_wkb(p1_position, POINT_DATA_SIZE, wkb_ndr, bin))
+      return 0;
+  }
+
+  bin->write_at_position(position, *poly_shapes);
+
+  return (uint) (opres - opres_orig);
+}
+
+
 uint Gis_polygon::init_from_wkb(const char *wkb, uint len, wkbByteOrder bo,
                                 String *res)
 {
@@ -746,9 +845,7 @@ uint Gis_polygon::init_from_wkb(const char *wkb, uint len, wkbByteOrder bo,
   if (len < 4)
     return 0;
 
-  if (!(n_linear_rings= wkb_get_uint(wkb, bo)))
-    return 0;
-
+  n_linear_rings= wkb_get_uint(wkb, bo);
   if (res->reserve(4, 512))
     return 0;
   wkb+= 4;
@@ -794,7 +891,7 @@ bool Gis_polygon::get_data_as_wkt(String *txt, const char **end) const
       return 1;
     n_points= uint4korr(data);
     data+= 4;
-    if (no_data(data, (SIZEOF_STORED_DOUBLE*2) * n_points) ||
+    if (no_data(data, POINT_DATA_SIZE * n_points) ||
 	txt->reserve(2 + ((MAX_DIGITS_IN_DOUBLE + 1) * 2 + 1) * n_points))
       return 1;
     txt->qs_append('(');
@@ -848,16 +945,16 @@ int Gis_polygon::area(double *ar, const char **end_of_data) const
     if (no_data(data, 4))
       return 1;
     n_points= uint4korr(data);
-    if (no_data(data, (SIZEOF_STORED_DOUBLE*2) * n_points))
+    if (no_data(data, POINT_DATA_SIZE * n_points))
       return 1;
     get_point(&prev_x, &prev_y, data+4);
-    data+= (4+SIZEOF_STORED_DOUBLE*2);
+    data+= (4+POINT_DATA_SIZE);
 
     while (--n_points)				// One point is already read
     {
       double x, y;
       get_point(&x, &y, data);
-      data+= (SIZEOF_STORED_DOUBLE*2);
+      data+= POINT_DATA_SIZE;
       lr_area+= (prev_x + x)* (prev_y - y);
       prev_x= x;
       prev_y= y;
@@ -884,7 +981,7 @@ int Gis_polygon::exterior_ring(String *result) const
   n_points= uint4korr(data);
   data+= 4;
   length= n_points * POINT_DATA_SIZE;
-  if (no_data(data, length) || result->reserve(1+4+4+ length))
+  if (no_data(data, length) || result->reserve(1 + 4 + 4 + length))
     return 1;
 
   result->q_append((char) wkb_ndr);
@@ -930,7 +1027,7 @@ int Gis_polygon::interior_ring_n(uint32 num, String *result) const
   n_points= uint4korr(data);
   points_size= n_points * POINT_DATA_SIZE;
   data+= 4;
-  if (no_data(data, points_size) || result->reserve(1+4+4+ points_size))
+  if (no_data(data, points_size) || result->reserve(1 + 4 + 4 + points_size))
     return 1;
 
   result->q_append((char) wkb_ndr);
@@ -969,16 +1066,16 @@ int Gis_polygon::centroid_xy(double *x, double *y) const
       return 1;
     org_n_points= n_points= uint4korr(data);
     data+= 4;
-    if (no_data(data, (SIZEOF_STORED_DOUBLE*2) * n_points))
+    if (no_data(data, POINT_DATA_SIZE * n_points))
       return 1;
     get_point(&prev_x, &prev_y, data);
-    data+= (SIZEOF_STORED_DOUBLE*2);
+    data+= POINT_DATA_SIZE;
 
     while (--n_points)				// One point is already read
     {
       double tmp_x, tmp_y;
       get_point(&tmp_x, &tmp_y, data);
-      data+= (SIZEOF_STORED_DOUBLE*2);
+      data+= POINT_DATA_SIZE;
       cur_area+= (prev_x + tmp_x) * (prev_y - tmp_y);
       cur_cx+= tmp_x;
       cur_cy+= tmp_y;
@@ -1018,6 +1115,49 @@ int Gis_polygon::centroid(String *result) const
   return create_point(result, x, y);
 }
 
+
+int Gis_polygon::store_shapes(Gcalc_shape_transporter *trn) const
+{
+  uint32 n_linear_rings;
+  const char *data= m_data;
+
+  if (trn->start_poly())
+    return 1;
+
+  if (no_data(data, 4))
+    return 1;
+  n_linear_rings= uint4korr(data);
+  data+= 4;
+
+  while (n_linear_rings--)
+  {
+    uint32 n_points;
+
+    if (no_data(data, 4))
+      return 1;
+    n_points= uint4korr(data);
+    data+= 4;
+    if (!n_points || no_data(data, POINT_DATA_SIZE * n_points))
+      return 1;
+
+    trn->start_ring();
+    while (--n_points)
+    {
+      double x, y;
+      get_point(&x, &y, data);
+      data+= POINT_DATA_SIZE;
+      if (trn->add_point(x, y))
+        return 1;
+    }
+    data+= POINT_DATA_SIZE;
+    trn->complete_ring();
+  }
+
+  trn->complete_poly();
+  return 0;
+}
+
+
 const Geometry::Class_info *Gis_polygon::get_class_info() const
 {
   return &polygon_class;
@@ -1046,7 +1186,7 @@ bool Gis_multi_point::init_from_wkt(Gis_read_stream *trs, String *wkb)
 
   for (;;)
   {
-    if (wkb->reserve(1+4, 512))
+    if (wkb->reserve(1 + 4, 512))
       return 1;
     wkb->q_append((char) wkb_ndr);
     wkb->q_append((uint32) wkb_point);
@@ -1061,6 +1201,32 @@ bool Gis_multi_point::init_from_wkt(Gis_read_stream *trs, String *wkb)
 }
 
 
+uint Gis_multi_point::init_from_opresult(String *bin,
+                                         const char *opres, uint32 n_shapes)
+{
+  uint bin_size, opres_size;
+  Gis_point p;
+  const char *opres_end;
+
+  bin_size= n_shapes * (WKB_HEADER_SIZE + POINT_DATA_SIZE) + 4;
+  opres_size= n_shapes * (4 + 8*2);
+ 
+  if (bin->reserve(bin_size, 512))
+    return 0;
+    
+  bin->q_append(n_shapes);
+  opres_end= opres + opres_size;
+  for (; opres < opres_end; opres+= (4 + 8*2))
+  {
+    bin->q_append((char)wkb_ndr);
+    bin->q_append((uint32)wkb_point);
+    if (!p.init_from_wkb(opres + 4, POINT_DATA_SIZE, wkb_ndr, bin))
+      return 0;
+  }
+  return opres_size;
+}
+
+
 uint Gis_multi_point::init_from_wkb(const char *wkb, uint len, wkbByteOrder bo,
                                     String *res)
 {
@@ -1099,7 +1265,7 @@ bool Gis_multi_point::get_data_as_wkt(String *txt, const char **end) const
 
   n_points= uint4korr(m_data);
   if (no_data(m_data+4,
-	      n_points * (SIZEOF_STORED_DOUBLE * 2 + WKB_HEADER_SIZE)) ||
+	      n_points * (POINT_DATA_SIZE + WKB_HEADER_SIZE)) ||
       txt->reserve(((MAX_DIGITS_IN_DOUBLE + 1) * 2 + 1) * n_points))
     return 1;
   *end= append_points(txt, n_points, m_data+4, WKB_HEADER_SIZE);
@@ -1140,6 +1306,35 @@ int Gis_multi_point::geometry_n(uint32 num, String *result) const
   return 0;
 }
 
+
+int Gis_multi_point::store_shapes(Gcalc_shape_transporter *trn) const
+{
+  uint32 n_points;
+  Gis_point pt;
+  const char *data= m_data;
+
+  if (no_data(data, 4))
+    return 1;
+  n_points= uint4korr(data);
+  data+= 4;
+
+  if (trn->start_collection(n_points))
+    return 1;
+
+  while (n_points--)
+  {
+    if (no_data(data, WKB_HEADER_SIZE))
+      return 1;
+    data+= WKB_HEADER_SIZE;
+    pt.set_data_ptr(data, (uint32) (m_data_end - data));
+    if (pt.store_shapes(trn))
+      return 1;
+    data+= pt.get_data_size();
+  }
+  return 0;
+}
+
+
 const Geometry::Class_info *Gis_multi_point::get_class_info() const
 {
   return &multipoint_class;
@@ -1182,10 +1377,9 @@ bool Gis_multi_line_string::init_from_wkt(Gis_read_stream *trs, String *wkb)
   {
     Gis_line_string ls;
 
-    if (wkb->reserve(1+4, 512))
+    if (wkb->reserve(1 + 4, 512))
       return 1;
-    wkb->q_append((char) wkb_ndr);
-    wkb->q_append((uint32) wkb_linestring);
+    wkb->q_append((char) wkb_ndr); wkb->q_append((uint32) wkb_linestring);
 
     if (trs->check_next_symbol('(') ||
 	ls.init_from_wkt(trs, wkb) ||
@@ -1200,15 +1394,44 @@ bool Gis_multi_line_string::init_from_wkt(Gis_read_stream *trs, String *wkb)
 }
 
 
+uint Gis_multi_line_string::init_from_opresult(String *bin,
+                                               const char *opres,
+                                               uint32 n_shapes)
+{
+  const char *opres_orig= opres;
+
+  if (bin->reserve(4, 512))
+    return 0;
+  bin->q_append(n_shapes);
+  
+  while (n_shapes--)
+  {
+    Gis_line_string ls;
+    int ls_len;
+
+    if (bin->reserve(WKB_HEADER_SIZE, 512))
+      return 0;
+
+    bin->q_append((char) wkb_ndr);
+    bin->q_append((uint32) wkb_linestring);
+
+    if (!(ls_len= ls.init_from_opresult(bin, opres)))
+      return 0;
+    opres+= ls_len;
+  }
+  return (uint) (opres - opres_orig);
+}
+
+
 uint Gis_multi_line_string::init_from_wkb(const char *wkb, uint len,
                                           wkbByteOrder bo, String *res)
 {
   uint32 n_line_strings;
   const char *wkb_orig= wkb;
 
-  if (len < 4)
+  if (len < 4 ||
+      (n_line_strings= wkb_get_uint(wkb, bo))< 1)
     return 0;
-  n_line_strings= wkb_get_uint(wkb, bo);
 
   if (res->reserve(4, 512))
     return 0;
@@ -1256,7 +1479,7 @@ bool Gis_multi_line_string::get_data_as_wkt(String *txt,
       return 1;
     n_points= uint4korr(data + WKB_HEADER_SIZE);
     data+= WKB_HEADER_SIZE + 4;
-    if (no_data(data, n_points * (SIZEOF_STORED_DOUBLE*2)) ||
+    if (no_data(data, n_points * POINT_DATA_SIZE) ||
 	txt->reserve(2 + ((MAX_DIGITS_IN_DOUBLE + 1) * 2 + 1) * n_points))
       return 1;
     txt->qs_append('(');
@@ -1386,6 +1609,35 @@ int Gis_multi_line_string::is_closed(int *closed) const
   return 0;
 }
 
+
+int Gis_multi_line_string::store_shapes(Gcalc_shape_transporter *trn) const
+{
+  uint32 n_lines;
+  Gis_line_string ls;
+  const char *data= m_data;
+
+  if (no_data(data, 4))
+    return 1;
+  n_lines= uint4korr(data);
+  data+= 4;
+
+  if (trn->start_collection(n_lines))
+    return 1;
+
+  while (n_lines--)
+  {
+    if (no_data(data, WKB_HEADER_SIZE))
+      return 1;
+    data+= WKB_HEADER_SIZE;
+    ls.set_data_ptr(data, (uint32) (m_data_end - data));
+    if (ls.store_shapes(trn))
+      return 1;
+    data+= ls.get_data_size();
+  }
+  return 0;
+}
+
+
 const Geometry::Class_info *Gis_multi_line_string::get_class_info() const
 {
   return &multilinestring_class;
@@ -1436,7 +1688,7 @@ bool Gis_multi_polygon::init_from_wkt(Gis_read_stream *trs, String *wkb)
 
   for (;;)  
   {
-    if (wkb->reserve(1+4, 512))
+    if (wkb->reserve(1 + 4, 512))
       return 1;
     wkb->q_append((char) wkb_ndr);
     wkb->q_append((uint32) wkb_polygon);
@@ -1491,6 +1743,37 @@ uint Gis_multi_polygon::init_from_wkb(const char *wkb, uint len,
 }
 
 
+uint Gis_multi_polygon::init_from_opresult(String *bin,
+                                           const char *opres, uint32 n_shapes)
+{
+  Gis_polygon p;
+  const char *opres_orig= opres;
+  uint p_len;
+  uint poly_shapes;
+  uint n_poly= 0;
+  uint32 np_pos= bin->length();
+
+  if (bin->reserve(4, 512))
+    return 0;
+    
+  bin->q_append(n_shapes);
+  while (n_shapes)
+  {
+    if (bin->reserve(1 + 4, 512))
+      return 0;
+    bin->q_append((char)wkb_ndr);
+    bin->q_append((uint32)wkb_polygon);
+    if (!(p_len= p.priv_init_from_opresult(bin, opres, n_shapes, &poly_shapes)))
+      return 0;
+    n_shapes-= poly_shapes;
+    opres+= p_len;
+    n_poly++;
+  }
+  bin->write_at_position(np_pos, n_poly);
+  return opres - opres_orig;
+}
+
+
 bool Gis_multi_polygon::get_data_as_wkt(String *txt, const char **end) const
 {
   uint32 n_polygons;
@@ -1517,7 +1800,7 @@ bool Gis_multi_polygon::get_data_as_wkt(String *txt, const char **end) const
         return 1;
       uint32 n_points= uint4korr(data);
       data+= 4;
-      if (no_data(data, (SIZEOF_STORED_DOUBLE * 2) * n_points) ||
+      if (no_data(data, POINT_DATA_SIZE * n_points) ||
 	  txt->reserve(2 + ((MAX_DIGITS_IN_DOUBLE + 1) * 2 + 1) * n_points,
 		       512))
 	return 1;
@@ -1678,6 +1961,35 @@ int Gis_multi_polygon::centroid(String *result) const
   return create_point(result, res_cx, res_cy);
 }
 
+
+int Gis_multi_polygon::store_shapes(Gcalc_shape_transporter *trn) const
+{
+  uint32 n_polygons;
+  Gis_polygon p;
+  const char *data= m_data;
+
+  if (no_data(data, 4))
+    return 1;
+  n_polygons= uint4korr(data);
+  data+= 4;
+
+  if (trn->start_collection(n_polygons))
+    return 1;
+
+  while (n_polygons--)
+  {
+    if (no_data(data, WKB_HEADER_SIZE))
+      return 1;
+    data+= WKB_HEADER_SIZE;
+    p.set_data_ptr(data, (uint32) (m_data_end - data));
+    if (p.store_shapes(trn))
+      return 1;
+    data+= p.get_data_size();
+  }
+  return 0;
+}
+
+
 const Geometry::Class_info *Gis_multi_polygon::get_class_info() const
 {
   return &multipolygon_class;
@@ -1749,6 +2061,45 @@ bool Gis_geometry_collection::init_from_wkt(Gis_read_stream *trs, String *wkb)
 }
 
 
+uint Gis_geometry_collection::init_from_opresult(String *bin,
+                                                 const char *opres,
+                                                 uint32 n_shapes)
+{
+  const char *opres_orig= opres;
+  Geometry_buffer buffer;
+  Geometry *geom;
+  int g_len;
+  uint32 wkb_type;
+
+  if (bin->reserve(4, 512))
+    return 0;
+  bin->q_append(n_shapes);
+  
+  while (n_shapes--)
+  {
+    switch ((Gcalc_function::shape_type) uint4korr(opres))
+    {
+      case Gcalc_function::shape_point:   wkb_type= wkb_point; break;
+      case Gcalc_function::shape_line:    wkb_type= wkb_linestring; break;
+      case Gcalc_function::shape_polygon: wkb_type= wkb_polygon; break;
+      default: wkb_type= 0; DBUG_ASSERT(FALSE);
+    };
+
+    if (bin->reserve(WKB_HEADER_SIZE, 512))
+      return 0;
+
+    bin->q_append((char) wkb_ndr);
+    bin->q_append(wkb_type);
+
+    if (!(geom= create_by_typeid(&buffer, wkb_type)) ||
+        !(g_len= geom->init_from_opresult(bin, opres, 1)))
+      return 0;
+    opres+= g_len;
+  }
+  return (uint) (opres - opres_orig);
+}
+
+
 uint Gis_geometry_collection::init_from_wkb(const char *wkb, uint len,
                                             wkbByteOrder bo, String *res)
 {
@@ -1898,7 +2249,7 @@ int Gis_geometry_collection::geometry_n(uint32 num, String *result) const
   } while (--num);
 
   /* Copy found object to result */
-  if (result->reserve(1+4+length))
+  if (result->reserve(1 + 4 + length))
     return 1;
   result->q_append((char) wkb_ndr);
   result->q_append((uint32) wkb_type);
@@ -1959,6 +2310,42 @@ bool Gis_geometry_collection::dimension(uint32 *res_dim, const char **end) const
   return 0;
 }
 
+
+int Gis_geometry_collection::store_shapes(Gcalc_shape_transporter *trn) const
+{
+  uint32 n_objects;
+  const char *data= m_data;
+  Geometry_buffer buffer;
+  Geometry *geom;
+
+  if (no_data(data, 4))
+    return 1;
+  n_objects= uint4korr(data);
+  data+= 4;
+
+  if (trn->start_collection(n_objects))
+    return 1;
+
+  while (n_objects--)
+  {
+    uint32 wkb_type;
+
+    if (no_data(data, WKB_HEADER_SIZE))
+      return 1;
+    wkb_type= uint4korr(data + 1);
+    data+= WKB_HEADER_SIZE;
+    if (!(geom= create_by_typeid(&buffer, wkb_type)))
+      return 1;
+    geom->set_data_ptr(data, (uint32) (m_data_end - data));
+    if (geom->store_shapes(trn))
+      return 1;
+
+    data+= geom->get_data_size();
+  }
+  return 0;
+}
+
+
 const Geometry::Class_info *Gis_geometry_collection::get_class_info() const
 {
   return &geometrycollection_class;
diff --git a/sql/spatial.h b/sql/spatial.h
index f778acd6c34..6f6c39b5b42 100644
--- a/sql/spatial.h
+++ b/sql/spatial.h
@@ -1,4 +1,4 @@
-/* Copyright (C) 2002-2006 MySQL AB
+/* Copyright (c) 2002, 2010, Oracle and/or its affiliates. All rights reserved.
 
    This program is free software; you can redistribute it and/or modify
    it under the terms of the GNU General Public License as published by
@@ -10,14 +10,19 @@
    GNU General Public License for more details.
 
    You should have received a copy of the GNU General Public License
-   along with this program; if not, write to the Free Software
-   Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA */
+   along with this program; if not, write to the Free Software Foundation,
+   51 Franklin Street, Suite 500, Boston, MA 02110-1335 USA */
 
 #ifndef _spatial_h
 #define _spatial_h
 
+#include "sql_string.h"                         /* String, LEX_STRING */
+#include <my_compiler.h>
+
 #ifdef HAVE_SPATIAL
 
+#include "gcalc_tools.h"
+
 const uint SRID_SIZE= 4;
 const uint SIZEOF_STORED_DOUBLE= 8;
 const uint POINT_DATA_SIZE= SIZEOF_STORED_DOUBLE*2; 
@@ -242,10 +247,13 @@ public:
   virtual const Class_info *get_class_info() const=0;
   virtual uint32 get_data_size() const=0;
   virtual bool init_from_wkt(Gis_read_stream *trs, String *wkb)=0;
-
   /* returns the length of the wkb that was read */
   virtual uint init_from_wkb(const char *wkb, uint len, wkbByteOrder bo,
                              String *res)=0;
+  virtual uint init_from_opresult(String *bin,
+                                  const char *opres, uint32 n_shapes=1)
+  { return init_from_wkb(opres + 4, UINT_MAX32, wkb_ndr, bin) + 4; }
+
   virtual bool get_data_as_wkt(String *txt, const char **end) const=0;
   virtual bool get_mbr(MBR *mbr, const char **end) const=0;
   virtual bool dimension(uint32 *dim, const char **end) const=0;
@@ -264,16 +272,20 @@ public:
   virtual int point_n(uint32 num, String *result) const { return -1; }
   virtual int interior_ring_n(uint32 num, String *result) const { return -1; }
   virtual int geometry_n(uint32 num, String *result) const { return -1; }
+  virtual int store_shapes(Gcalc_shape_transporter *trn) const=0;
 
 public:
   static Geometry *create_by_typeid(Geometry_buffer *buffer, int type_id);
+
   static Geometry *construct(Geometry_buffer *buffer,
                              const char *data, uint32 data_len);
   static Geometry *create_from_wkt(Geometry_buffer *buffer,
 				   Gis_read_stream *trs, String *wkt,
 				   bool init_stream=1);
-  static Geometry *create_from_wkb(Geometry_buffer *buffer, const char *wkb, 
-                                   uint32 len, String *res);
+  static Geometry *create_from_wkb(Geometry_buffer *buffer,
+                                   const char *wkb, uint32 len, String *res);
+  static int create_from_opresult(Geometry_buffer *g_buf,
+                                  String *res, Gcalc_result_receiver &rr);
   int as_wkt(String *wkt, const char **end)
   {
     uint32 len= (uint) get_class_info()->m_name.length;
@@ -369,6 +381,7 @@ public:
     *end= 0;					/* No default end */
     return 0;
   }
+  int store_shapes(Gcalc_shape_transporter *trn) const;
   const Class_info *get_class_info() const;
 };
 
@@ -397,6 +410,7 @@ public:
     *end= 0;					/* No default end */
     return 0;
   }
+  int store_shapes(Gcalc_shape_transporter *trn) const;
   const Class_info *get_class_info() const;
 };
 
@@ -411,6 +425,13 @@ public:
   uint32 get_data_size() const;
   bool init_from_wkt(Gis_read_stream *trs, String *wkb);
   uint init_from_wkb(const char *wkb, uint len, wkbByteOrder bo, String *res);
+  uint priv_init_from_opresult(String *bin, const char *opres,
+                               uint32 n_shapes, uint32 *poly_shapes);
+  uint init_from_opresult(String *bin, const char *opres, uint32 n_shapes)
+  {
+    uint32 foo;
+    return priv_init_from_opresult(bin, opres, n_shapes, &foo);
+  }
   bool get_data_as_wkt(String *txt, const char **end) const;
   bool get_mbr(MBR *mbr, const char **end) const;
   int area(double *ar, const char **end) const;
@@ -425,6 +446,7 @@ public:
     *end= 0;					/* No default end */
     return 0;
   }
+  int store_shapes(Gcalc_shape_transporter *trn) const;
   const Class_info *get_class_info() const;
 };
 
@@ -439,6 +461,7 @@ public:
   uint32 get_data_size() const;
   bool init_from_wkt(Gis_read_stream *trs, String *wkb);
   uint init_from_wkb(const char *wkb, uint len, wkbByteOrder bo, String *res);
+  uint init_from_opresult(String *bin, const char *opres, uint32 n_shapes);
   bool get_data_as_wkt(String *txt, const char **end) const;
   bool get_mbr(MBR *mbr, const char **end) const;
   int num_geometries(uint32 *num) const;
@@ -449,6 +472,7 @@ public:
     *end= 0;					/* No default end */
     return 0;
   }
+  int store_shapes(Gcalc_shape_transporter *trn) const;
   const Class_info *get_class_info() const;
 };
 
@@ -463,6 +487,7 @@ public:
   uint32 get_data_size() const;
   bool init_from_wkt(Gis_read_stream *trs, String *wkb);
   uint init_from_wkb(const char *wkb, uint len, wkbByteOrder bo, String *res);
+  uint init_from_opresult(String *bin, const char *opres, uint32 n_shapes);
   bool get_data_as_wkt(String *txt, const char **end) const;
   bool get_mbr(MBR *mbr, const char **end) const;
   int num_geometries(uint32 *num) const;
@@ -475,6 +500,7 @@ public:
     *end= 0;					/* No default end */
     return 0;
   }
+  int store_shapes(Gcalc_shape_transporter *trn) const;
   const Class_info *get_class_info() const;
 };
 
@@ -501,7 +527,9 @@ public:
     *end= 0;					/* No default end */
     return 0;
   }
+  int store_shapes(Gcalc_shape_transporter *trn) const;
   const Class_info *get_class_info() const;
+  uint init_from_opresult(String *bin, const char *opres, uint32 n_shapes);
 };
 
 
@@ -515,11 +543,13 @@ public:
   uint32 get_data_size() const;
   bool init_from_wkt(Gis_read_stream *trs, String *wkb);
   uint init_from_wkb(const char *wkb, uint len, wkbByteOrder bo, String *res);
+  uint init_from_opresult(String *bin, const char *opres, uint32 n_shapes);
   bool get_data_as_wkt(String *txt, const char **end) const;
   bool get_mbr(MBR *mbr, const char **end) const;
   int num_geometries(uint32 *num) const;
   int geometry_n(uint32 num, String *result) const;
   bool dimension(uint32 *dim, const char **end) const;
+  int store_shapes(Gcalc_shape_transporter *trn) const;
   const Class_info *get_class_info() const;
 };