diff --git a/Common/Context.h b/Common/Context.h
index 0dabad2cb5e87410f192aa66d79746175a2bb6a3..8d22c2ebd9e0f630a9915d4401149199ae3bf178 100644
--- a/Common/Context.h
+++ b/Common/Context.h
@@ -183,6 +183,7 @@ public :
   struct{
     int format;
     int eps_quality, eps_background, eps_font_size;
+    int eps_occlusion_culling, eps_best_root;
     char *eps_font;
     double eps_line_width_factor, eps_point_size_factor;
     int jpeg_quality;
diff --git a/Common/DefaultOptions.h b/Common/DefaultOptions.h
index 1066e3cf92142546565261084c317347272a2e84..3b6f3dfe45f162199a8201ce0ef0876b2b4db957 100644
--- a/Common/DefaultOptions.h
+++ b/Common/DefaultOptions.h
@@ -975,6 +975,10 @@ StringXNumber PrintOptions_Number[] = {
     "Size factor for points in postscript output" },
   { F|O, "EpsQuality" , opt_print_eps_quality , 1 ,
     "Postscript quality (1=simple sort, 2=recursive sort)" },
+  { F|O, "EpsOcclusionCulling" , opt_print_eps_occlusion_culling , 0 ,
+    "Cull occluded primitives (to reduce PostScript file size)" },
+  { F|O, "EpsBestRoot" , opt_print_eps_best_root , 0 ,
+    "Try to minimize primitive splitting in recursive sort" },
 
   { F|O, "Format" , opt_print_format , FORMAT_AUTO , 
     "File format (10=automatic)" }, 
diff --git a/Common/Options.cpp b/Common/Options.cpp
index 692c97c8e3175c667f9c79cc28511e3e38a46551..45562a7f4d5796de74d3ffa5b91ac23a98894573 100644
--- a/Common/Options.cpp
+++ b/Common/Options.cpp
@@ -1,4 +1,4 @@
-// $Id: Options.cpp,v 1.96 2002-11-17 17:08:32 geuzaine Exp $
+// $Id: Options.cpp,v 1.97 2002-12-11 17:37:17 geuzaine Exp $
 //
 // Copyright (C) 1997 - 2002 C. Geuzaine, J.-F. Remacle
 //
@@ -3009,6 +3009,16 @@ double opt_print_eps_quality(OPT_ARGS_NUM){
     CTX.print.eps_quality = (int)val;
   return CTX.print.eps_quality;
 }
+double opt_print_eps_occlusion_culling(OPT_ARGS_NUM){
+  if(action & GMSH_SET) 
+    CTX.print.eps_occlusion_culling = (int)val;
+  return CTX.print.eps_occlusion_culling;
+}
+double opt_print_eps_best_root(OPT_ARGS_NUM){
+  if(action & GMSH_SET) 
+    CTX.print.eps_best_root = (int)val;
+  return CTX.print.eps_best_root;
+}
 double opt_print_eps_background(OPT_ARGS_NUM){
   if(action & GMSH_SET) 
     CTX.print.eps_background = (int)val;
diff --git a/Common/Options.h b/Common/Options.h
index f5dd32b3fbfb0212151dc16378eefaf8b0183880..2581c9f4921f9577a257ace4778a7e4b5262e162 100644
--- a/Common/Options.h
+++ b/Common/Options.h
@@ -455,6 +455,8 @@ double opt_view_point_type(OPT_ARGS_NUM);
 double opt_view_line_type(OPT_ARGS_NUM);
 double opt_print_format(OPT_ARGS_NUM);
 double opt_print_eps_quality(OPT_ARGS_NUM);
+double opt_print_eps_occlusion_culling(OPT_ARGS_NUM);
+double opt_print_eps_best_root(OPT_ARGS_NUM);
 double opt_print_eps_background(OPT_ARGS_NUM);
 double opt_print_eps_font_size(OPT_ARGS_NUM);
 double opt_print_eps_line_width_factor(OPT_ARGS_NUM);
diff --git a/Graphics/CreateFile.cpp b/Graphics/CreateFile.cpp
index 88580e9241440d9a623ec23736405c9ff7f88087..6f33d836121d22c97ebe2621f6f9174bc2289cd5 100644
--- a/Graphics/CreateFile.cpp
+++ b/Graphics/CreateFile.cpp
@@ -1,4 +1,4 @@
-// $Id: CreateFile.cpp,v 1.35 2002-10-11 07:14:34 geuzaine Exp $
+// $Id: CreateFile.cpp,v 1.36 2002-12-11 17:37:17 geuzaine Exp $
 //
 // Copyright (C) 1997 - 2002 C. Geuzaine, J.-F. Remacle
 //
@@ -184,6 +184,8 @@ void CreateOutputFile (char *name, int format) {
     psformat = (format==FORMAT_PS || format==FORMAT_PSTEX) ? GL2PS_PS : GL2PS_EPS;
     pssort = (CTX.print.eps_quality==1) ? GL2PS_SIMPLE_SORT : GL2PS_BSP_SORT;
     psoptions = GL2PS_SIMPLE_LINE_OFFSET | GL2PS_SILENT |
+      (CTX.print.eps_occlusion_culling ? GL2PS_OCCLUSION_CULL : 0) |
+      (CTX.print.eps_best_root ? GL2PS_BEST_ROOT : 0) |
       (CTX.print.eps_background ? GL2PS_DRAW_BACKGROUND : 0) |
       (format==FORMAT_PSTEX ? GL2PS_NO_TEXT : 0) |
       (format==FORMAT_EPSTEX ? GL2PS_NO_TEXT : 0);
diff --git a/Graphics/gl2ps.cpp b/Graphics/gl2ps.cpp
index 1cf12afcc1dbfae773bc817b186ed3ca4121fce4..6d0927449275354d682ebe6a592b7aee6203a2bf 100644
--- a/Graphics/gl2ps.cpp
+++ b/Graphics/gl2ps.cpp
@@ -2,7 +2,7 @@
  * GL2PS, an OpenGL to PostScript Printing Library
  * Copyright (C) 1999-2002  Christophe Geuzaine 
  *
- * $Id: gl2ps.cpp,v 1.52 2002-11-17 02:30:13 geuzaine Exp $
+ * $Id: gl2ps.cpp,v 1.53 2002-12-11 17:37:17 geuzaine Exp $
  *
  * E-mail: geuz@geuz.org
  * URL: http://www.geuz.org/gl2ps/
@@ -33,7 +33,7 @@
 /* The gl2ps context. gl2ps is not thread safe (we should create a
    local GL2PScontext during gl2psBeginPage). */
 
-GL2PScontext *gl2ps=NULL;
+GL2PScontext *gl2ps = NULL;
 
 /* Some 'system' utility routines */
 
@@ -128,23 +128,23 @@ void *gl2psListPointer(GL2PSlist *list, GLint index){
 }
 
 void gl2psListSort(GL2PSlist *list,
-		     int (*fcmp)(const void *a, const void *b)){
+		   int (*fcmp)(const void *a, const void *b)){
   qsort(list->array, list->n, list->size, fcmp);
 }
 
 void gl2psListAction(GL2PSlist *list, 
-		       void (*action)(void *data, void *dummy)){
+		     void (*action)(void *data, void *dummy)){
   GLint i, dummy;
 
-  for(i=0 ; i<gl2psListNbr(list) ; i++)
+  for(i = 0; i < gl2psListNbr(list); i++)
     (*action)(gl2psListPointer(list, i), &dummy);
 }
 
 void gl2psListActionInverse(GL2PSlist *list, 
-			      void (*action)(void *data, void *dummy)){
+			    void (*action)(void *data, void *dummy)){
   GLint i, dummy;
 
-  for(i=gl2psListNbr(list) ; i>0 ; i--)
+  for(i = gl2psListNbr(list); i > 0; i--)
     (*action)(gl2psListPointer(list, i-1), &dummy);
 }
 
@@ -179,12 +179,14 @@ void gl2psGetNormal(GLfloat *a, GLfloat *b, GLfloat *c){
     c[1] = c[1] / norm;
     c[2] = c[2] / norm;
   }
-  else
-    gl2psMsg(GL2PS_WARNING, "Bad plane in BSP tree");
+  else{
+    c[0] = c[1] = 0.;
+    c[2] = 1.;
+  }
 }
 
 void gl2psGetPlane(GL2PSprimitive *prim, GL2PSplane plane){
-  GL2PSxyz v={0., 0., 0.}, w={0., 0., 0.};
+  GL2PSxyz v = {0., 0., 0.}, w = {0., 0., 0.};
 
   switch(prim->type){
   case GL2PS_TRIANGLE :
@@ -241,7 +243,7 @@ void gl2psGetPlane(GL2PSprimitive *prim, GL2PSplane plane){
 }
 
 void gl2psCutEdge(GL2PSvertex *a, GL2PSvertex *b, GL2PSplane plane, 
-		    GL2PSvertex *c){
+		  GL2PSvertex *c){
   GL2PSxyz v;
   GLfloat  sect;
 
@@ -261,8 +263,8 @@ void gl2psCutEdge(GL2PSvertex *a, GL2PSvertex *b, GL2PSplane plane,
 }
 
 void gl2psCreateSplittedPrimitive(GL2PSprimitive *parent, GL2PSplane plane,
-				    GL2PSprimitive *child, GLshort numverts,
-				    GLshort *index0, GLshort *index1){
+				  GL2PSprimitive *child, GLshort numverts,
+				  GLshort *index0, GLshort *index1){
   GLshort i;
 
   if(numverts > 4){
@@ -282,7 +284,7 @@ void gl2psCreateSplittedPrimitive(GL2PSprimitive *parent, GL2PSplane plane,
   child->numverts = numverts;
   child->verts = (GL2PSvertex *)gl2psMalloc(numverts * sizeof(GL2PSvertex));
 
-  for(i=0 ; i<numverts ; i++){
+  for(i = 0; i < numverts; i++){
     if(index1[i] < 0)
       child->verts[i] = parent->verts[index0[i]];
     else
@@ -292,10 +294,10 @@ void gl2psCreateSplittedPrimitive(GL2PSprimitive *parent, GL2PSplane plane,
 }
 
 void gl2psAddIndex(GLshort *index0, GLshort *index1, GLshort *nb, 
-		     GLshort i, GLshort j){
+		   GLshort i, GLshort j){
   GLint k;
 
-  for(k=0 ; k<*nb ; k++)
+  for(k = 0; k < *nb; k++)
     if((index0[k] == i && index1[k] == j) ||
        (index1[k] == i && index0[k] == j)) return;
 
@@ -309,18 +311,18 @@ GLshort gl2psGetIndex(GLshort i, GLshort num){
 }
 
 GLint gl2psTestSplitPrimitive(GL2PSprimitive *prim, GL2PSplane plane){
-  GLint     type=GL2PS_COINCIDENT;
+  GLint     type = GL2PS_COINCIDENT;
   GLshort   i, j;
   GLfloat   d[5]; 
 
-  for(i = 0 ; i < prim->numverts ; i++){  
+  for(i = 0; i < prim->numverts; i++){  
     d[i] = gl2psComparePointPlane(prim->verts[i].xyz, plane);
   }
 
   if(prim->type == GL2PS_POINT)
     return 0;
   else{
-    for(i = 0 ; i < prim->numverts ; i++){
+    for(i = 0; i < prim->numverts; i++){
       j = gl2psGetIndex(i, prim->numverts);
       if(d[j] > GL2PS_EPSILON){
 	if(type == GL2PS_COINCIDENT)      type = GL2PS_IN_BACK_OF;
@@ -345,7 +347,7 @@ GLint gl2psSplitPrimitive(GL2PSprimitive *prim, GL2PSplane plane,
 
   type = GL2PS_COINCIDENT;
 
-  for(i = 0 ; i < prim->numverts ; i++){  
+  for(i = 0; i < prim->numverts; i++){  
     d[i] = gl2psComparePointPlane(prim->verts[i].xyz, plane);
   }
 
@@ -356,7 +358,7 @@ GLint gl2psSplitPrimitive(GL2PSprimitive *prim, GL2PSplane plane,
     else                           type = GL2PS_COINCIDENT;
     break;
   default :
-    for(i = 0 ; i < prim->numverts ; i++){
+    for(i = 0; i < prim->numverts; i++){
       j = gl2psGetIndex(i, prim->numverts);
       if(d[j] > GL2PS_EPSILON){
 	if(type == GL2PS_COINCIDENT)      type = GL2PS_IN_BACK_OF;
@@ -397,7 +399,7 @@ GLint gl2psSplitPrimitive(GL2PSprimitive *prim, GL2PSplane plane,
 }
 
 void gl2psDivideQuad(GL2PSprimitive *quad, 
-		       GL2PSprimitive **t1, GL2PSprimitive **t2){
+		     GL2PSprimitive **t1, GL2PSprimitive **t2){
   *t1 = (GL2PSprimitive*)gl2psMalloc(sizeof(GL2PSprimitive));
   *t2 = (GL2PSprimitive*)gl2psMalloc(sizeof(GL2PSprimitive));
   (*t1)->type = (*t2)->type = GL2PS_TRIANGLE;
@@ -441,17 +443,22 @@ int gl2psTrianglesFirst(const void *a, const void *b){
 }
 
 GLint gl2psFindRoot(GL2PSlist *primitives, GL2PSprimitive **root){
-  GLint          i, j, count, best=1000000, index=0;
+  GLint          i, j, count, best = 1000000, index = 0;
   GL2PSprimitive *prim1, *prim2;
   GL2PSplane     plane;
+  GLint maxp;
 
   if(gl2ps->options & GL2PS_BEST_ROOT){
     *root = *(GL2PSprimitive**)gl2psListPointer(primitives, 0);
-    for(i=0 ; i<gl2psListNbr(primitives) ; i++){
+    maxp = gl2psListNbr(primitives);
+    if(maxp > gl2ps->maxbestroot){
+      maxp = gl2ps->maxbestroot;
+    }
+    for(i = 0; i < maxp; i++){
       prim1 = *(GL2PSprimitive**)gl2psListPointer(primitives, i);
       gl2psGetPlane(prim1, plane);
-      count=0;
-      for(j=0 ; j<gl2psListNbr(primitives) ; j++){
+      count = 0;
+      for(j = 0; j < gl2psListNbr(primitives); j++){
 	if(j != i){
 	  prim2 = *(GL2PSprimitive**)gl2psListPointer(primitives, j);
 	  count += gl2psTestSplitPrimitive(prim2, plane); 
@@ -475,7 +482,7 @@ GLint gl2psFindRoot(GL2PSlist *primitives, GL2PSprimitive **root){
 }
 
 void gl2psFreePrimitive(void *a, void *b){
-  GL2PSprimitive *q ;
+  GL2PSprimitive *q;
   
   q = *(GL2PSprimitive**)a;
   gl2psFree(q->verts);
@@ -542,10 +549,10 @@ void gl2psBuildBspTree(GL2PSbsptree *tree, GL2PSlist *primitives){
   frontlist = gl2psListCreate(1, 2, sizeof(GL2PSprimitive*));
   backlist = gl2psListCreate(1, 2, sizeof(GL2PSprimitive*));
 
-  for(i=0 ; i<gl2psListNbr(primitives) ; i++){
+  for(i = 0; i < gl2psListNbr(primitives); i++){
     if(i != index){
       prim = *(GL2PSprimitive**)gl2psListPointer(primitives,i);
-      switch(gl2psSplitPrimitive(prim,tree->plane,&frontprim,&backprim)){
+      switch(gl2psSplitPrimitive(prim, tree->plane, &frontprim, &backprim)){
       case GL2PS_COINCIDENT:
 	gl2psAddPrimitiveInList(prim, tree->primitives);
 	break;
@@ -586,9 +593,9 @@ void gl2psBuildBspTree(GL2PSbsptree *tree, GL2PSlist *primitives){
   gl2psListDelete(primitives);
 }
 
-void  gl2psTraverseBspTree(GL2PSbsptree *tree, GL2PSxyz eye, GLfloat epsilon,
-			     GLboolean (*compare)(GLfloat f1, GLfloat f2),
-			     void (*action)(void *data, void *dummy)){
+void gl2psTraverseBspTree(GL2PSbsptree *tree, GL2PSxyz eye, GLfloat epsilon,
+			  GLboolean (*compare)(GLfloat f1, GLfloat f2),
+			  void (*action)(void *data, void *dummy)){
   GLfloat result;
 
   if(!tree) return;
@@ -603,17 +610,302 @@ void  gl2psTraverseBspTree(GL2PSbsptree *tree, GL2PSxyz eye, GLfloat epsilon,
   else if(compare(-epsilon, result)){ 
     gl2psTraverseBspTree(tree->front, eye, epsilon, compare, action);
     gl2psListAction(tree->primitives, action);
-    gl2psTraverseBspTree(tree->back, eye, epsilon, compare, action) ;
+    gl2psTraverseBspTree(tree->back, eye, epsilon, compare, action);
   }
   else{
     gl2psTraverseBspTree(tree->front, eye, epsilon, compare, action);
-    gl2psTraverseBspTree(tree->back, eye, epsilon, compare, action) ;
+    gl2psTraverseBspTree(tree->back, eye, epsilon, compare, action);
   }
 }
 
-/* Boundary contruction */
+/* The 2D sorting routines (for occlusion culling) */
+
+GLint gl2psGetPlaneFromPoints(GL2PSxyz a, GL2PSxyz b, GL2PSplane plane){  
+  GLfloat  n; 
+  plane[0] = b[1] - a[1];
+  plane[1] = a[0] - b[0];
+  n = sqrt(plane[0]*plane[0] + plane[1]*plane[1]);
+  plane[2]=0.;
+  if(n != 0.){
+    plane[0] /= n;
+    plane[1] /= n;
+    plane[3] = -plane[0]*a[0]-plane[1]*a[1]; 
+    return 1;
+  }
+  else{
+    plane[0] = -1.0;
+    plane[1] = 0.;
+    plane[3] = a[0];
+    return 0;
+  }
+}
 
-#define GL2PS_BOUNDARY_OFFSET 0
+void gl2psFreeBspImageTree(GL2PSbsptree2d **tree){
+  if(*tree){
+    if((*tree)->back)  gl2psFreeBspImageTree(&(*tree)->back);
+    if((*tree)->front) gl2psFreeBspImageTree(&(*tree)->front);
+    gl2psFree(*tree);
+    *tree = NULL;
+  }
+}
+
+GLint gl2psCheckPoint(GL2PSxyz point, GL2PSplane plane){
+  GLfloat pt_dis;
+  pt_dis = gl2psComparePointPlane(point, plane);
+  if(pt_dis > GL2PS_EPSILON)        return GL2PS_POINT_INFRONT;
+  else if(pt_dis < -GL2PS_EPSILON)  return GL2PS_POINT_BACK;
+  else                              return GL2PS_POINT_COINCIDENT;
+}
+
+void gl2psAddPlanesInBspTreeImage(GL2PSprimitive *prim,
+				  GL2PSbsptree2d **tree){
+  GLint ret = 0;
+  GLint i;
+  GLint offset = 0;
+  GL2PSbsptree2d *head = NULL, *cur = NULL;
+
+  if((*tree == NULL) && (prim->numverts > 2)){
+    head=(GL2PSbsptree2d*)gl2psMalloc(sizeof(GL2PSbsptree2d));
+    for(i = 0; i < prim->numverts-1; i++){
+      if(!gl2psGetPlaneFromPoints(prim->verts[i].xyz,
+                                  prim->verts[i+1].xyz,
+				  head->plane)){
+        if(prim->numverts-i > 3)  
+	  offset++;
+	else{
+	  gl2psFree(head);
+	  return;
+	}
+      }
+      else break;
+    }
+    head->back = NULL;
+    head->front = NULL;
+    for(i = 2+offset; i < prim->numverts; i++){
+      ret = gl2psCheckPoint(prim->verts[i].xyz, head->plane);
+      if(ret != 0) break;
+    }
+    switch(ret){
+    case GL2PS_POINT_INFRONT :
+      cur = head;
+      for(i = 1+offset; i < (prim->numverts-1); i++){
+        if(cur->front == NULL)
+          cur->front = (GL2PSbsptree2d*)gl2psMalloc(sizeof(GL2PSbsptree2d));
+        if(gl2psGetPlaneFromPoints(prim->verts[i].xyz,
+				   prim->verts[i+1].xyz,
+				   cur->front->plane)){
+	  cur = cur->front;
+	  cur->front = NULL;
+	  cur->back = NULL;
+	}
+      }
+      if(cur->front == NULL)
+        cur->front = (GL2PSbsptree2d*)gl2psMalloc(sizeof(GL2PSbsptree2d));
+      if(gl2psGetPlaneFromPoints(prim->verts[i].xyz,
+                                 prim->verts[offset].xyz,
+				 cur->front->plane)){
+	cur->front->front = NULL;
+	cur->front->back = NULL;
+      }
+      else{
+        gl2psFree(cur->front);
+	cur = NULL;
+      }
+      break;
+    case GL2PS_POINT_BACK :
+      for(i = 0; i < 4; i++)
+        head->plane[i] = -head->plane[i];
+      cur = head;
+      for(i = 1+offset; i < (prim->numverts-1); i++){
+        if(cur->front == NULL)
+          cur->front = (GL2PSbsptree2d*)gl2psMalloc(sizeof(GL2PSbsptree2d));
+        if(gl2psGetPlaneFromPoints(prim->verts[i+1].xyz,
+				   prim->verts[i].xyz,
+				   cur->front->plane)){
+	  cur = cur->front;
+	  cur->front = NULL;
+	  cur->back = NULL;
+	}
+      }
+      if(cur->front == NULL)
+        cur->front = (GL2PSbsptree2d*)gl2psMalloc(sizeof(GL2PSbsptree2d));
+      if(gl2psGetPlaneFromPoints(prim->verts[offset].xyz,
+                                 prim->verts[i].xyz,
+				 cur->front->plane)){
+	cur->front->front = NULL;
+	cur->front->back = NULL;
+      }
+      else{
+        gl2psFree(cur->front);
+	cur = NULL;
+      }
+      break;
+    default:
+      gl2psFree(head);
+      return;
+    }
+    (*tree) = head;
+  }
+}
+
+GLint gl2psCheckPrimitive(GL2PSprimitive *prim, GL2PSplane plane){
+  GLint i;
+  GLint pos;
+  pos = gl2psCheckPoint(prim->verts[0].xyz, plane);
+  for(i = 1; i < prim->numverts; i++){
+    pos |= gl2psCheckPoint(prim->verts[i].xyz, plane);
+    if(pos == (GL2PS_POINT_INFRONT | GL2PS_POINT_BACK)) return GL2PS_SPANNING;
+  }
+  if(pos & GL2PS_POINT_INFRONT)   return GL2PS_IN_FRONT_OF;
+  else if(pos & GL2PS_POINT_BACK) return GL2PS_IN_BACK_OF;
+  else                            return GL2PS_COINCIDENT;
+}
+
+GL2PSprimitive* gl2psCreateSplitPrimitive2D(GL2PSprimitive *parent,
+                                            GLint numverts,
+					    GL2PSvertex *vertx){
+  GLint i;
+  GL2PSprimitive *child = (GL2PSprimitive*)gl2psMalloc(sizeof(GL2PSprimitive));
+  switch(numverts){
+  case 1 : child->type = GL2PS_POINT; break;
+  case 2 : child->type = GL2PS_LINE; break;
+  case 3 : child->type = GL2PS_TRIANGLE; break;
+  case 4 : child->type = GL2PS_QUADRANGLE; break;
+  }
+  child->boundary = 0; /* not done! */
+  child->dash = parent->dash;
+  child->width = parent->width;
+  child->numverts = numverts;
+  child->verts = (GL2PSvertex *)gl2psMalloc(numverts * sizeof(GL2PSvertex));
+  for(i = 0; i < numverts; i++)
+    child->verts[i] = vertx[i];
+  return child;
+}
+
+void gl2psSplitPrimitive2D(GL2PSprimitive *prim, 
+			   GL2PSplane plane, 
+			   GL2PSprimitive **front, 
+			   GL2PSprimitive **back){
+
+  // cur will hold the position of current vertex
+  // prev will holds the position of previous vertex
+  // prev0 will holds the position of vertex number 0
+  // v1 and v2 represent the current and previous vertexs respectively
+  // flag will represents that should the current be checked against the plane
+  GLint cur = -1, prev = -1, i, v1 = 0, v2 = 0, flag = 1, prev0 = -1;
+  
+  // list of vertexs which will go in front and back Primitive
+  GL2PSvertex *front_list = NULL, *back_list = NULL;
+  
+  // number of vertex in front and back list
+  GLint front_count = 0, back_count = 0;
+
+  for(i = 0; i <= prim->numverts; i++){
+    v1 = i;
+    if(v1 == prim->numverts){
+      if(prim->numverts < 3) break;
+      v1 = 0;
+      v2 = prim->numverts-1;
+      cur = prev0;
+    }
+    else if(flag){
+      cur = gl2psCheckPoint(prim->verts[v1].xyz, plane);
+      if(i == 0)
+        prev0 = cur;
+    } 
+    if(((prev == -1) || (prev == cur) || (prev == 0) || (cur == 0)) &&
+       (i < prim->numverts)){
+      if(cur == GL2PS_POINT_INFRONT){
+	front_count++;
+	front_list = (GL2PSvertex*)gl2psRealloc(front_list,
+						sizeof(GL2PSvertex)*front_count);
+        front_list[front_count-1] = prim->verts[v1];
+      }
+      else if(cur == GL2PS_POINT_BACK){
+	back_count++;
+	back_list = (GL2PSvertex*)gl2psRealloc(back_list,
+					       sizeof(GL2PSvertex)*back_count);
+        back_list[back_count-1] = prim->verts[v1];
+      }
+      else{
+	front_count++;
+	front_list = (GL2PSvertex*)gl2psRealloc(front_list,
+						sizeof(GL2PSvertex)*front_count);
+        front_list[front_count-1] = prim->verts[v1];
+	back_count++;
+	back_list = (GL2PSvertex*)gl2psRealloc(back_list,
+					       sizeof(GL2PSvertex)*back_count);
+        back_list[back_count-1] = prim->verts[v1];
+      }
+      flag = 1;
+    }
+    else if((prev != cur) && (cur != 0) && (prev != 0)){
+      if(v1 != 0){
+	v2 = v1-1;
+	i--;
+      }
+      front_count++;
+      front_list = (GL2PSvertex*)gl2psRealloc(front_list,
+					      sizeof(GL2PSvertex)*front_count);
+      gl2psCutEdge(&prim->verts[v2],
+                   &prim->verts[v1],
+		   plane,
+		   &front_list[front_count-1]);
+      back_count++;
+      back_list = (GL2PSvertex*)gl2psRealloc(back_list,
+					     sizeof(GL2PSvertex)*back_count);
+      back_list[back_count-1] = front_list[front_count-1];
+      flag = 0;
+    }
+    prev = cur;
+  }
+  *front = gl2psCreateSplitPrimitive2D(prim, front_count, front_list);
+  *back = gl2psCreateSplitPrimitive2D(prim, back_count, back_list);
+  gl2psFree(front_list);
+  gl2psFree(back_list);
+}
+
+GLint gl2psAddInImageTree(GL2PSprimitive *prim, GL2PSbsptree2d **tree){
+  GLint ret = 0;
+  GL2PSprimitive *frontprim = NULL, *backprim = NULL;
+  if(*tree == NULL){
+    gl2psAddPlanesInBspTreeImage(prim, tree);
+    return 1;
+  }
+  else{
+    switch(gl2psCheckPrimitive(prim, (*tree)->plane)){
+    case GL2PS_IN_BACK_OF: return gl2psAddInImageTree(prim, &(*tree)->back);
+    case GL2PS_IN_FRONT_OF: 
+      if((*tree)->front != NULL) return gl2psAddInImageTree(prim, &(*tree)->front);
+      else                       return 0;
+    case GL2PS_SPANNING:
+      gl2psSplitPrimitive2D(prim, (*tree)->plane, &frontprim, &backprim);
+      ret = gl2psAddInImageTree(backprim, &(*tree)->back);
+      if((*tree)->front != NULL)
+        if(gl2psAddInImageTree(frontprim, &(*tree)->front))  
+	  ret = 1;
+      gl2psFree(frontprim->verts);
+      gl2psFree(frontprim);
+      gl2psFree(backprim->verts);
+      gl2psFree(backprim);
+      return ret;
+    case GL2PS_COINCIDENT:
+      if(prim->numverts < 3) return 1;
+      else                   return 0;
+    }
+  }
+  return 0;
+}
+
+void gl2psAddInImage(void *a, void *b){
+  GL2PSprimitive *prim;
+  prim = *(GL2PSprimitive **)a;
+  if(gl2psAddInImageTree(prim, &gl2ps->image)){
+    prim->depth = -1.;
+  }
+}
+
+/* Boundary contruction */
 
 void gl2psAddBoundaryInList(GL2PSprimitive *prim, GL2PSlist *list){
   GL2PSprimitive *b;
@@ -621,14 +913,14 @@ void gl2psAddBoundaryInList(GL2PSprimitive *prim, GL2PSlist *list){
   GL2PSxyz        c;
 
   c[0] = c[1] = c[2] = 0.;
-  for(i=0 ; i<prim->numverts ; i++){
+  for(i = 0; i < prim->numverts; i++){
     c[0] += prim->verts[i].xyz[0];
     c[1] += prim->verts[i].xyz[1];
   }
   c[0] /= prim->numverts;
   c[1] /= prim->numverts;
 
-  for(i=0 ; i<prim->numverts ; i++){
+  for(i = 0; i < prim->numverts; i++){
     if(prim->boundary & (GLint)pow(2., i)){
       b = (GL2PSprimitive*)gl2psMalloc(sizeof(GL2PSprimitive));
       b->type = GL2PS_LINE;
@@ -638,6 +930,7 @@ void gl2psAddBoundaryInList(GL2PSprimitive *prim, GL2PSlist *list){
       b->numverts = 2;
       b->verts = (GL2PSvertex *)gl2psMalloc(2 * sizeof(GL2PSvertex));
 
+#define GL2PS_BOUNDARY_OFFSET 0
 #if GL2PS_BOUNDARY_OFFSET
       v[0] = c[0] - prim->verts[i].xyz[0];
       v[1] = c[1] - prim->verts[i].xyz[1];
@@ -679,14 +972,14 @@ void gl2psAddBoundaryInList(GL2PSprimitive *prim, GL2PSlist *list){
 
 }
 
-void  gl2psBuildPolygonBoundary(GL2PSbsptree *tree){
+void gl2psBuildPolygonBoundary(GL2PSbsptree *tree){
   GLint          i, n;
   GL2PSprimitive *prim;
 
   if(!tree) return;
   gl2psBuildPolygonBoundary(tree->back);
   n = gl2psListNbr(tree->primitives);
-  for(i=0 ; i<n ; i++){
+  for(i = 0; i < n; i++){
     prim = *(GL2PSprimitive**)gl2psListPointer(tree->primitives, i);
     if(prim->boundary) gl2psAddBoundaryInList(prim, tree->primitives);
   }
@@ -696,9 +989,9 @@ void  gl2psBuildPolygonBoundary(GL2PSbsptree *tree){
 /* The feedback buffer parser */
 
 void gl2psAddPolyPrimitive(GLshort type, GLshort numverts, 
-			     GL2PSvertex *verts, GLint offset, 
-			     GLshort dash, GLfloat width,
-			     GLshort boundary){
+			   GL2PSvertex *verts, GLint offset, 
+			   GLshort dash, GLfloat width,
+			   GLshort boundary){
   GLshort         i;
   GLfloat         factor, units, area, dZ, dZdX, dZdY, maxdZ;
   GL2PSprimitive *prim;
@@ -726,7 +1019,8 @@ void gl2psAddPolyPrimitive(GLshort type, GLshort numverts,
   }
   else if(offset && type == GL2PS_TRIANGLE){
 
-    /* needs some more work... */
+    /* This needs some more work... */
+
     if(gl2ps->sort == GL2PS_SIMPLE_SORT){    
       factor = gl2ps->offset[0];
       units = gl2ps->offset[1];
@@ -756,8 +1050,6 @@ void gl2psAddPolyPrimitive(GLshort type, GLshort numverts,
 
     dZ = factor * maxdZ + units;
 
-    /* printf("dZ = %g  (fact=%g  units=%g)\n", dZ, factor, units); */
-
     prim->verts[0].xyz[2] += dZ;
     prim->verts[1].xyz[2] += dZ;
     prim->verts[2].xyz[2] += dZ;
@@ -801,9 +1093,9 @@ GLint gl2psGetVertex(GL2PSvertex *v, GLfloat *p){
 }
 
 GLint gl2psParseFeedbackBuffer(void){
-  GLint        i, used, count, v, vtot, offset=0;
-  GLshort      boundary, flag, dash=0;
-  GLfloat      lwidth=1., psize=1.;
+  GLint        i, used, count, v, vtot, offset = 0;
+  GLshort      boundary, flag, dash = 0;
+  GLfloat      lwidth = 1., psize = 1.;
   GLfloat     *current;
   GL2PSvertex  vertices[3];
 
@@ -860,9 +1152,9 @@ GLint gl2psParseFeedbackBuffer(void){
 	vtot++;
 	if(v == 2){
 	  if(boundary){
-	    if(!count && vtot==2) flag = 1|2|4;
+	    if(!count && vtot == 2) flag = 1|2|4;
 	    else if(!count) flag = 2|4;
-	    else if(vtot==2) flag = 1|2;
+	    else if(vtot == 2) flag = 1|2;
 	    else flag = 2;
 	  }
 	  else
@@ -886,21 +1178,21 @@ GLint gl2psParseFeedbackBuffer(void){
       break;      
     case GL_PASS_THROUGH_TOKEN :
       switch((GLint)current[1]){
-      case GL2PS_BEGIN_POLYGON_OFFSET_FILL : offset=1; break;
-      case GL2PS_END_POLYGON_OFFSET_FILL : offset=0; break;
-      case GL2PS_BEGIN_POLYGON_BOUNDARY : boundary=1; break;
-      case GL2PS_END_POLYGON_BOUNDARY : boundary=0; break;
-      case GL2PS_BEGIN_LINE_STIPPLE : dash=4; break;
-      case GL2PS_END_LINE_STIPPLE : dash=0; break;
+      case GL2PS_BEGIN_POLYGON_OFFSET_FILL : offset = 1; break;
+      case GL2PS_END_POLYGON_OFFSET_FILL : offset = 0; break;
+      case GL2PS_BEGIN_POLYGON_BOUNDARY : boundary = 1; break;
+      case GL2PS_END_POLYGON_BOUNDARY : boundary = 0; break;
+      case GL2PS_BEGIN_LINE_STIPPLE : dash = 4; break;
+      case GL2PS_END_LINE_STIPPLE : dash = 0; break;
       case GL2PS_SET_POINT_SIZE : 
-	current+=2; 
-	used-=2; 
-	psize=current[1];
+	current += 2; 
+	used -= 2; 
+	psize = current[1];
 	break;
       case GL2PS_SET_LINE_WIDTH : 
-	current+=2; 
-	used-=2; 
-	lwidth=current[1];
+	current += 2; 
+	used -= 2; 
+	lwidth = current[1];
 	break;
       }
       current += 2; 
@@ -925,7 +1217,7 @@ GLboolean gl2psSameColor(GL2PSrgba rgba1, GL2PSrgba rgba2){
   
 GLboolean gl2psVertsSameColor(const GL2PSprimitive *prim){
   int i;
-  for(i=1; i<prim->numverts; i++)
+  for(i = 1; i < prim->numverts; i++)
     if(!gl2psSameColor(prim->verts[0].rgba, prim->verts[i].rgba))
       return 0;
   return 1;
@@ -1120,7 +1412,7 @@ void gl2psPrintPostScriptHeader(void){
 	  "%%%%Page: 1 1\n"
 	  "%%%%BeginPageSetup\n");
 
-  if (gl2ps->options & GL2PS_LANDSCAPE)
+  if(gl2ps->options & GL2PS_LANDSCAPE)
     fprintf(gl2ps->stream,
 	    "%d 0 translate 90 rotate\n",
 	    viewport[3]);
@@ -1249,7 +1541,7 @@ void gl2psPrintPostScriptFooter(void){
 	  "%%%%EOF\n");
 }
 
-/* The LaTeX routines. */
+/* The LaTeX routines */
 
 void gl2psPrintTeXHeader(void){
   GLint   viewport[4];
@@ -1258,15 +1550,15 @@ void gl2psPrintTeXHeader(void){
 
   glGetIntegerv(GL_VIEWPORT, viewport);
 
-  if(gl2ps->filename && strlen(gl2ps->filename)<256){
-    for(i=strlen(gl2ps->filename)-1 ; i>=0 ; i--){
+  if(gl2ps->filename && strlen(gl2ps->filename) < 256){
+    for(i = strlen(gl2ps->filename)-1; i >= 0; i--){
       if(gl2ps->filename[i] == '.'){
 	strncpy(name, gl2ps->filename, i);
-	name[i]='\0';
+	name[i] = '\0';
 	break;
       }
     }
-    if(i<=0) strcpy(name, gl2ps->filename);
+    if(i <= 0) strcpy(name, gl2ps->filename);
   }
   else
     strcpy(name, "untitled");
@@ -1304,12 +1596,13 @@ void gl2psPrintTeXFooter(void){
 
 /* The public routines */
 
-GL2PSDLL_API void gl2psBeginPage(char *title, char *producer, 
-				   GLint format, GLint sort, GLint options, 
-				   GLint colormode, GLint colorsize, 
-				   GL2PSrgba *colormap, GLint buffersize, 
-				   FILE *stream, char *filename){
+GL2PSDLL_API void gl2psBeginPage(const char *title, const char *producer, 
+				 GLint format, GLint sort, GLint options, 
+				 GLint colormode, GLint colorsize, 
+				 GL2PSrgba *colormap, GLint buffersize, 
+				 FILE *stream, const char *filename){
   gl2ps = (GL2PScontext*)gl2psMalloc(sizeof(GL2PScontext));
+  gl2ps->maxbestroot = 10;
   gl2ps->format = format;
   gl2ps->title = title;
   gl2ps->producer = producer;
@@ -1327,6 +1620,7 @@ GL2PSDLL_API void gl2psBeginPage(char *title, char *producer,
   gl2ps->lastrgba[2] = -1.;
   gl2ps->lastrgba[3] = -1.;
   gl2ps->lastlinewidth = -1.;
+  gl2ps->image = NULL;
   gl2ps->primitives = gl2psListCreate(500, 500, sizeof(GL2PSprimitive*));
 
   if(gl2ps->colormode == GL_RGBA){
@@ -1354,7 +1648,7 @@ GL2PSDLL_API void gl2psBeginPage(char *title, char *producer,
 
 GL2PSDLL_API GLint gl2psEndPage(void){
   GL2PSbsptree   *root;
-  GL2PSxyz        eye={0., 0., 100000.};
+  GL2PSxyz        eye = {0., 0., 100000.};
   GLint           shademodel, res;
   void          (*phead)(void) = 0;
   void          (*pprim)(void *a, void *b) = 0;
@@ -1402,6 +1696,10 @@ GL2PSDLL_API GLint gl2psEndPage(void){
       break;
     case GL2PS_SIMPLE_SORT :
       gl2psListSort(gl2ps->primitives, gl2psCompareDepth);
+      if(gl2ps->options & GL2PS_OCCLUSION_CULL){
+	gl2psListAction(gl2ps->primitives,gl2psAddInImage);
+        gl2psFreeBspImageTree(&gl2ps->image);
+      }
       gl2psListActionInverse(gl2ps->primitives, pprim);
       gl2psListAction(gl2ps->primitives, gl2psFreePrimitive);
       gl2psListDelete(gl2ps->primitives);
@@ -1411,12 +1709,11 @@ GL2PSDLL_API GLint gl2psEndPage(void){
       root = (GL2PSbsptree*)gl2psMalloc(sizeof(GL2PSbsptree));
       gl2psBuildBspTree(root, gl2ps->primitives);
       if(gl2ps->boundary) gl2psBuildPolygonBoundary(root);
-      /* Occlusion culling is not implemented yet...
       if(gl2ps->options & GL2PS_OCCLUSION_CULL){
 	gl2psTraverseBspTree(root, eye, -(float)GL2PS_EPSILON, gl2psLess,
 			     gl2psAddInImage);
+	gl2psFreeBspImageTree(&gl2ps->image);
       }
-      */
       gl2psTraverseBspTree(root, eye, (float)GL2PS_EPSILON, gl2psGreater, 
 			   pprim);
       gl2psFreeBspTree(root);
@@ -1427,8 +1724,9 @@ GL2PSDLL_API GLint gl2psEndPage(void){
     }
     pfoot();
     fflush(gl2ps->stream);
-
+    
   }
+  
 
   if(gl2ps->colormap) gl2psFree(gl2ps->colormap);
 
@@ -1438,7 +1736,7 @@ GL2PSDLL_API GLint gl2psEndPage(void){
   return res;
 }
 
-GL2PSDLL_API void gl2psText(char *str, char *fontname, GLint fontsize){
+GL2PSDLL_API void gl2psText(const char *str, const char *fontname, GLint fontsize){
   GLfloat         pos[4];
   GL2PSprimitive  *prim;
   GLboolean       valid;
diff --git a/Graphics/gl2ps.h b/Graphics/gl2ps.h
index e6510f4cc949a1e3749f8259810655309769c3ed..cfb5ea30531915c4924c23bad23811a4285eaf71 100644
--- a/Graphics/gl2ps.h
+++ b/Graphics/gl2ps.h
@@ -2,7 +2,7 @@
  * GL2PS, an OpenGL to PostScript Printing Library
  * Copyright (C) 1999-2002  Christophe Geuzaine
  *
- * $Id: gl2ps.h,v 1.29 2002-11-12 19:11:50 geuzaine Exp $
+ * $Id: gl2ps.h,v 1.30 2002-12-11 17:37:17 geuzaine Exp $
  *
  * E-mail: geuz@geuz.org
  * URL: http://www.geuz.org/gl2ps/
@@ -55,7 +55,7 @@
 #endif /* __APPLE__ */
 
 
-#define GL2PS_VERSION                    0.63
+#define GL2PS_VERSION                    0.7
 #define GL2PS_NONE                       0
 
 /* Output file format */
@@ -123,6 +123,12 @@
 #define GL2PS_IN_BACK_OF                 3
 #define GL2PS_SPANNING                   4
 
+/* 2D BSP tree primitive comparison */
+
+#define GL2PS_POINT_COINCIDENT           0
+#define GL2PS_POINT_INFRONT              1
+#define GL2PS_POINT_BACK                 2
+
 /* Pass through options */
 
 #define GL2PS_BEGIN_POLYGON_OFFSET_FILL  1
@@ -138,6 +144,13 @@ typedef GLfloat GL2PSrgba[4];
 typedef GLfloat GL2PSxyz[3];
 typedef GLfloat GL2PSplane[4];
 
+typedef struct _GL2PSbsptree2d GL2PSbsptree2d;
+
+struct _GL2PSbsptree2d {
+  GL2PSplane plane;
+  GL2PSbsptree2d *front, *back;
+};
+
 typedef struct {
   GLint nmax, size, incr, n;
   char *array;
@@ -169,13 +182,14 @@ typedef struct {
 } GL2PSprimitive;
 
 typedef struct {
-  GLint format, sort, options, colorsize, colormode, buffersize;
-  char *title, *producer, *filename;
+  GLint format, sort, options, colorsize, colormode, buffersize, maxbestroot;
+  const char *title, *producer, *filename;
   GLboolean shade, boundary;
   GLfloat *feedback, offset[2];
   GL2PSrgba *colormap, lastrgba, threshold;
   float lastlinewidth;
   GL2PSlist *primitives;
+  GL2PSbsptree2d *image;
   FILE *stream;
 } GL2PScontext;
 
@@ -186,13 +200,13 @@ typedef struct {
 extern "C" {
 #endif
 
-GL2PSDLL_API void  gl2psBeginPage(char *title, char *producer, 
+GL2PSDLL_API void  gl2psBeginPage(const char *title, const char *producer, 
 				  GLint format, GLint sort, GLint options, 
 				  GLint colormode, GLint colorsize, 
 				  GL2PSrgba *colormap, GLint buffersize, 
-				  FILE *stream, char *filename);
+				  FILE *stream, const char *filename);
 GL2PSDLL_API GLint gl2psEndPage(void);
-GL2PSDLL_API void  gl2psText(char *str, char *fontname, GLint size);
+GL2PSDLL_API void  gl2psText(const char *str, const char *fontname, GLint size);
 GL2PSDLL_API void  gl2psEnable(GLint mode);
 GL2PSDLL_API void  gl2psDisable(GLint mode);
 GL2PSDLL_API void  gl2psPointSize(GLfloat value);