1 #include "gp.h"
2 #include "induct.h"
3 
path_through_nonuni_gp(nodein,nodeout,plane)4 SPATH *path_through_nonuni_gp(nodein, nodeout, plane)
5      NODES *nodein, *nodeout;
6      GROUNDPLANE *plane;
7 {
8   G_nodes *in, *out;
9   Nonuni_gp *gp = plane->nonuni;
10   SPATH *path, *get_a_nonuni_path();
11 
12   if (nodein->gp_node != NULL)
13     in = nodein->gp_node;
14   else
15     GP_PANIC("path_through_nonuni_gp: nodein->gp_node == NULL!");
16 
17   if (nodeout->gp_node != NULL)
18     out = nodeout->gp_node;
19   else
20     GP_PANIC("path_through_nonuni_gp: nodeout->gp_node == NULL!");
21 
22   if (in == out) {
23     GP_PANIC("path_through_nonuni_gp: Path starts and ends at same node!");
24     /*
25     fprintf(stderr, "Warning: user gp node at (%lg, %lg, %lg) (meters) and node at (%lg, %lg, %lg)\n", nodein->x, nodein->y, nodein->z, nodeout->x,
26 	    nodeout->y, nodeout->z);
27     fprintf(stderr, "in Nonuniform gp %s correspond to the same node\n",
28 	    gp->grndp->name);
29     return NULL;
30     */
31   }
32 
33   clear_nonuni_marks(gp->nodelist);
34   path = get_a_nonuni_path(in, out, gp, NULL);
35 
36   /* if z-directed segs, add them to list */
37   if (gp->num_z_pts != 1)
38     GP_PANIC("find_non_uni_gp: You haven't added z-segs to path");
39 
40   return path;
41 }
42 
43 /* clear all the flags as NOT_CHECKED */
clear_nonuni_marks(node)44 clear_nonuni_marks(node)
45      G_nodes *node;
46 {
47   while(node != NULL) {
48     node->flag = NOT_CHECKED;
49     node = node->next;
50   }
51 }
52 
53 /*find the nearest node in the leaf cell that contains this global xyz point */
find_nearest_nonuni_node(xg,yg,zg,gp)54 G_nodes *find_nearest_nonuni_node(xg, yg, zg, gp)
55      double xg, yg, zg;
56      Nonuni_gp *gp;
57 {
58   Gcell *nearest;
59   G_nodes *node;
60 
61   double x,y,z;
62 
63   get_nonuni_coords(xg, yg, zg, gp, &x, &y, &z);
64 
65   /* find the cell this point is in */
66   nearest = get_containing_cell(x, y, gp->root_cell);
67 
68   if (!is_in_cell(x,y,nearest))
69     return NULL;
70 
71   /* now find nearest node along the edges of the cell */
72   return find_nearest_edge_node(x, y, nearest);
73 
74 }
75 
get_a_nonuni_path(node,node_goal,plane,nodes_so_far)76 SPATH *get_a_nonuni_path(node, node_goal, plane, nodes_so_far)
77      G_nodes *node, *node_goal;
78      Nonuni_gp *plane;
79      Llist *nodes_so_far;
80 {
81 #define MAXchoices 4
82   nonuni_choice_list choices[MAXchoices];
83 
84   int nchoices;
85   SPATH *path;
86   seg_ptr sptr;
87   double distance_to_goal;
88   double new_dist;
89   G_nodes *neighbor;
90   int i;
91 
92   distance_to_goal = get_node_dist(node_goal->x, node_goal->y, node);
93 
94   nchoices = 0;
95 
96   sptr.type = NORMAL;
97 
98   if (node->adjacent[N] != NULL && node->n_segs != NULL) {
99     neighbor = node->adjacent[N];
100     new_dist = get_node_dist(node_goal->x, node_goal->y, neighbor);
101     nchoices += add_nonuni_choice(&choices[nchoices], nodes_so_far,
102 		      node->n_segs, neighbor,
103 		      distance_to_goal - new_dist);
104   }
105   if (node->adjacent[E] != NULL && node->e_segs != NULL) {
106     neighbor = node->adjacent[E];
107     new_dist = get_node_dist(node_goal->x, node_goal->y, neighbor);
108     nchoices += add_nonuni_choice(&choices[nchoices], nodes_so_far,
109 		      node->e_segs, neighbor,
110 		      distance_to_goal - new_dist);
111   }
112   if (node->adjacent[S] != NULL && node->adjacent[S]->n_segs != NULL) {
113     neighbor = node->adjacent[S];
114     new_dist = get_node_dist(node_goal->x, node_goal->y, neighbor);
115     nchoices += add_nonuni_choice(&choices[nchoices], nodes_so_far,
116 		      node->adjacent[S]->n_segs, neighbor,
117 		      distance_to_goal - new_dist);
118   }
119   if (node->adjacent[W] != NULL && node->adjacent[W]->e_segs != NULL) {
120     neighbor = node->adjacent[W];
121     new_dist = get_node_dist(node_goal->x, node_goal->y, neighbor);
122     nchoices += add_nonuni_choice(&choices[nchoices], nodes_so_far,
123 		      node->adjacent[W]->e_segs, neighbor,
124 		      distance_to_goal - new_dist);
125   }
126 
127   /* sort in descreasing order by rank */
128   sort_nonuni_choices(choices, nchoices);
129 
130   /* add current node to list of passed through nodes */
131   nodes_so_far = add_ptr_to_list( (void *)node, nodes_so_far);
132 
133   for(path = NULL, i = 0; path == NULL && i < nchoices; i++) {
134     sptr.segp = (void *)choices[i].seg;
135     if (node_goal == choices[i].node) {
136       path = add_seg_to_list(sptr, path);
137       /* increment_usage(choices[i].seg); not used anymore (for overlap) */
138     }
139     else {
140       /* call myself with next node */
141       path = get_a_nonuni_path(choices[i].node, node_goal, plane,
142 			       nodes_so_far);
143       if (path != NULL) {
144 	/* there is a path, so lets add this seg to it */
145 	path = add_seg_to_list(sptr,path);
146 	/* increment_usage(choices[i].seg); not used anymore */
147       }
148     }
149   }
150 
151   /* free this node */
152   free(nodes_so_far);
153 
154   if (path == NULL)
155     /* mark it is a node from which no path was found */
156     node->flag = CHECKED;
157 
158   return path;
159 
160 
161 }
162 
163 /* sort choices DECREASING by rank.  Not the most efficient */
sort_nonuni_choices(choices,num)164 sort_nonuni_choices(choices, num)
165      nonuni_choice_list *choices;
166      int num;
167 {
168   int i, j;
169   nonuni_choice_list temp;
170 
171   for(i = 0; i < num - 1; i++)
172     for(j = num - 1; j > i; j--)
173       if (choices[j-1].rank < choices[j].rank) {
174 	temp = choices[j-1];
175 	choices[j-1] = choices[j];
176 	choices[j] = temp;
177       }
178 }
179 
180 
add_ptr_to_list(ptr,list)181 Llist *add_ptr_to_list(ptr, list)
182      void *ptr;
183      Llist *list;
184 {
185   Llist *new_llist;
186 
187   new_llist = (Llist *)gp_malloc(sizeof(Llist));
188 
189   new_llist->ptr = ptr;
190   new_llist->next = list;
191 
192   return new_llist;
193 }
194 
195 /* add or not a possible direction to search and rank it */
add_nonuni_choice(choice,nodes_so_far,segs,node,delta)196 int add_nonuni_choice(choice, nodes_so_far, segs, node, delta)
197      nonuni_choice_list *choice;
198      Llist *nodes_so_far;
199      SEGMENT **segs;
200      G_nodes *node;
201      double delta;
202 {
203   if (segs != NULL && !is_ptr_in_list((void *)node, nodes_so_far)
204       && node->flag == NOT_CHECKED) {
205     choice->seg = segs[0];
206     choice->node = node;
207     choice->rank = delta;
208 
209     return 1;
210   }
211   else
212     return 0;
213 }
214 
is_ptr_in_list(ptr,list)215 is_ptr_in_list(ptr, list)
216      void *ptr;
217      Llist *list;
218 {
219   while (list != NULL) {
220     if (list->ptr == ptr)
221       return TRUE;
222     list = list->next;
223   }
224 
225   return FALSE;
226 }
227 
228 /* find the cell that contains this point */
get_containing_cell(x,y,cell)229 Gcell *get_containing_cell(x, y, cell)
230      double x,y;
231      Gcell *cell;
232 {
233 
234   switch(get_children_type(cell)) {
235   case NONE:
236     /* we are a leaf and the point is in us */
237     return cell;
238   case BI:
239     return get_containing_bi_cell( x, y, cell);
240     break;
241   case GRID_2D:
242     return get_containing_grid_cell( x, y, cell );
243     break;
244   default:
245     GP_PANIC("Unknown child type in set_cell_coords")
246     break;
247   }
248 }
249 
250 /* figure out which of the grid of children contains (x,y) */
get_containing_grid_cell(x,y,cell)251 Gcell *get_containing_grid_cell( x, y, cell )
252      double x, y;
253      Gcell *cell;
254 {
255   double xn, yn;  /* normalized coords */
256   Grid_2d *grid;
257   int x_i, y_i;
258 
259   grid = (Grid_2d *)cell->children;
260 
261   /* you've written a function to do this: get_grid_indices() */
262 
263   xn = (x - cell->x0)/(cell->x1 - cell->x0);
264   yn = (y - cell->y0)/(cell->y1 - cell->y0);
265 
266   if (xn >= 1)
267     x_i = grid->x_cells - 1;
268   else if (xn < 0)
269     x_i = 0;
270   else
271     /* round to the nearest cell */
272     x_i = ((int) (xn * grid->x_cells));
273 
274   if (yn >= 1)
275     y_i = grid->y_cells - 1;
276   else if (yn < 0)
277     y_i = 0;
278   else
279     /* round to the nearest cell */
280     y_i = ((int) (yn * grid->y_cells));
281 
282   return get_containing_cell(x,y,grid->kids[grid->y_cells - y_i - 1][x_i]);
283 }
284 
285 
286 /* figure out which of the two children of cell contains (x,y) */
get_containing_bi_cell(x,y,cell)287 Gcell *get_containing_bi_cell( x, y, cell)
288      double x, y;
289      Gcell *cell;
290 {
291   double xn, yn;  /* normalized coords */
292   Bi *two_kids;
293 
294   two_kids = (Bi *)cell->children;
295 
296   xn = (x - cell->x0)/(cell->x1 - cell->x0);
297   yn = (y - cell->y0)/(cell->y1 - cell->y0);
298 
299   if (get_bi_type(two_kids) == NS) {
300     if (yn > 0.5)
301       return get_containing_cell(x,y,two_kids->child1);
302     else
303       return get_containing_cell(x,y,two_kids->child2);
304   }
305   else if (get_bi_type(two_kids) == EW) {
306     if (xn > 0.5)
307       return get_containing_cell(x,y,two_kids->child1);
308     else
309       return get_containing_cell(x,y,two_kids->child2);
310   }
311   else
312     GP_PANIC("get_containing bi_kids: Unknown bi_type!");
313 }
314 
is_in_cell(x,y,cell)315 is_in_cell(x, y, cell)
316      double x,y;
317      Gcell *cell;
318 {
319   if (x < cell->x0 || x > cell->x1 || y < cell->y0 || y > cell->y1)
320     return FALSE;
321   else
322     return TRUE;
323 }
324 
325 /* looks for the nearest node to x,y along the edges of the cell */
326 /* Called by find_nearest_nonuni_node() and make_equiv_rect() */
find_nearest_edge_node(x,y,cell)327 G_nodes *find_nearest_edge_node(x, y, cell)
328      double x,y;
329      Gcell *cell;
330 {
331   double xn, yn;
332   double x0 = get_x0(cell);
333   double y0 = get_y0(cell);
334   double x1 = get_x1(cell);
335   double y1 = get_y1(cell);
336   int is_y_bigger;
337 
338   G_nodes **nodes;
339   int quad;  /* quadrant */
340   char direction;
341 
342   if (!is_leaf(cell))
343     GP_PANIC("find_nearest_edge_node: This isn't a leaf cell!");
344 
345   nodes = &(cell->bndry.nodes[0]);
346 
347   /* put in normalized coords: 1x1 square with origin at center */
348 
349   xn = (x - x0) / (x1 - x0)  - 0.5;
350   yn = (y - y0) / (y1 - y0)  - 0.5;
351 
352   /* which node is this closest to */
353 
354   if (xn >= 0 && yn >= 0)
355     quad = NE;
356   else if (xn <= 0 && yn >= 0)
357     quad = NW;
358   else if (xn <= 0 && yn <= 0)
359     quad = SW;
360   else if (xn >= 0 && yn <= 0)
361     quad = SE;
362 
363   /* which direction points to the next nearest node */
364 
365   if (fabs(xn) < fabs(yn))
366     is_y_bigger = TRUE;
367   else
368     is_y_bigger = FALSE;
369 
370   /* scan along the given direction to find the nearest node */
371   if (quad == NE) {
372     if (is_y_bigger == TRUE)
373       direction = W;
374     else
375       direction = S;
376   }
377   else if (quad == NW) {
378     if (is_y_bigger == TRUE)
379       direction = E;
380     else
381       direction = S;
382   }
383   else if (quad == SW) {
384     if (is_y_bigger == TRUE)
385       direction = E;
386     else
387       direction = N;
388   }
389   else if (quad == SE) {
390     if (is_y_bigger == TRUE)
391       direction = W;
392     else
393       direction = N;
394   }
395 
396   return scan_edge(x, y, nodes[quad], direction);
397 }
398 
scan_edge(x,y,node,dir)399 G_nodes *scan_edge(x, y, node, dir)
400      double x,y;
401      G_nodes *node;
402      char dir;
403 {
404 
405   /* scan adjacency direction from node on until we minimize distance */
406 
407   double dist, tmp;
408 
409   dist = get_node_dist(x,y,node);
410 
411   /* there better always be an adjacent node since the adjacent
412      node to this one should be farther than the original node! */
413   if (node->adjacent[dir] == NULL)
414     GP_PANIC("scan_edge: adjacent[dir] == NULL!");
415 
416   while( (tmp = get_node_dist(x,y,node->adjacent[dir])) < dist) {
417     dist = tmp;
418     node = node->adjacent[dir];
419   }
420 
421   return node;
422 }
423 
get_node_dist(x,y,node)424 double get_node_dist(x,y,node)
425      double x,y;
426      G_nodes *node;
427 {
428   return get_dist(x - node->x, y - node->y);
429 }
430 
get_dist(x,y)431 double get_dist(x,y)
432      double x,y;
433 {
434   return sqrt(x*x + y*y);
435 }
436 
make_nonuni_Mlist(plane,pMlist)437 int make_nonuni_Mlist(plane, pMlist)
438      GROUNDPLANE *plane;
439      MELEMENT **pMlist;
440 {
441   int counter = 0;
442 
443   /* make meshes of the kids */
444   make_children_meshes(plane->nonuni->root_cell, pMlist, &counter);
445 
446   if (plane->nonuni->num_z_pts != 1) {
447     GP_PANIC("make_nonuni_Mlist: z_mesh function not written!");
448     /* make_z_meshes(plane->nonuni, pMlist, &counter); */
449   }
450 
451   if (counter > plane->numesh)
452     /* can be less (from edge holes), just not greater */
453     GP_PANIC("make_nonuni_Mlist: number of meshes too big!");
454 
455   return counter;
456 }
457 
make_children_meshes(cell,pMlist,pcount)458 make_children_meshes(cell, pMlist, pcount)
459      Gcell *cell;
460      MELEMENT **pMlist;
461      int *pcount;
462 {
463   switch (get_children_type(cell)) {
464   case NONE:
465     *pcount += make_leaf_mesh(cell, &(pMlist[*pcount]));
466     break;
467   case BI:
468     make_children_meshes( ( (Bi *)cell->children )->child1, pMlist, pcount);
469     make_children_meshes( ( (Bi *)cell->children )->child2, pMlist, pcount);
470     break;
471   case GRID_2D:
472     make_grid_children_meshes( (Grid_2d *)cell->children, pMlist, pcount);
473     break;
474   default:
475     GP_PANIC("make_children_meshes: Unknown child type!");
476     break;
477   }
478 }
479 
make_grid_children_meshes(grid,pMlist,pcount)480 make_grid_children_meshes( grid, pMlist, pcount)
481      Grid_2d *grid;
482      MELEMENT **pMlist;
483      int *pcount;
484 {
485   int i, j;
486 
487   for(i = 0; i < grid->y_cells; i++)
488     for(j = 0; j < grid->x_cells; j++)
489       make_children_meshes( grid->kids[i][j], pMlist, pcount);
490 }
491 
make_leaf_mesh(cell,pMlist)492 make_leaf_mesh(cell, pMlist)
493      Gcell *cell;
494      MELEMENT **pMlist;
495 {
496   G_nodes **nodes;
497   int bad;
498 
499   *pMlist = NULL;
500 
501   if (!is_leaf(cell))
502     GP_PANIC("make_leaf_mesh: not a leaf");
503 
504   nodes = cell->bndry.nodes;
505 
506   if (nodes[SW]->n_segs == NULL || nodes[SW]->e_segs == NULL
507       || nodes[NW]->e_segs == NULL || nodes[SE]->n_segs == NULL) {
508     /* this is a hole on an edge.  If it's not on an edge, we
509        must report an error since we haven't implemented a hole consisting
510        of multiple cells (only one mesh for the lot) */
511     if (!is_hole(cell))
512       GP_PANIC("make_leaf_mesh: NULL segments bordering non-hole cell!");
513 
514     bad = FALSE;
515     if (nodes[SW]->n_segs == NULL && nodes[SW]->cells[NE] != NULL
516         && nodes[SW]->cells[NW] != NULL)
517       bad = TRUE;
518     if (nodes[SW]->e_segs == NULL && nodes[SW]->cells[NE] != NULL
519         && nodes[SW]->cells[SE] != NULL)
520       bad = TRUE;
521     if (nodes[NW]->e_segs == NULL && nodes[NW]->cells[NE] != NULL
522         && nodes[NW]->cells[SE] != NULL)
523       bad = TRUE;
524     if (nodes[SE]->n_segs == NULL && nodes[SE]->cells[NE] != NULL
525         && nodes[SE]->cells[NW] != NULL)
526       bad = TRUE;
527 
528     if (bad == TRUE)
529       GP_PANIC("make_leaf_mesh: Adjacent hole cells not allowed!");
530 
531     /* if not bad, then do nothing */
532     return 0;
533   }
534   else {
535     /* make a mesh as required */
536     *pMlist = add_edge_segs_to_list(nodes[SW], nodes[NW], N, 1, *pMlist);
537     *pMlist = add_edge_segs_to_list(nodes[NW], nodes[NE], E, 1, *pMlist);
538     *pMlist = add_edge_segs_to_list(nodes[SE], nodes[NE], N, -1, *pMlist);
539     *pMlist = add_edge_segs_to_list(nodes[SW], nodes[SE], E, -1, *pMlist);
540 
541     return 1;
542   }
543 
544 }
545 
add_edge_segs_to_list(node0,node1,dir,signofelem,mfirst)546 MELEMENT *add_edge_segs_to_list(node0, node1, dir, signofelem, mfirst)
547      G_nodes *node0, *node1;
548      char dir;
549      int signofelem;
550      MELEMENT *mfirst;
551 {
552   G_nodes *node;
553   MELEMENT *plist = mfirst;
554   FILAMENT *fil;
555   MELEMENT *melem;
556 
557   for(node = node0; node != node1; node = node->adjacent[dir]) {
558     if (dir == N)
559       fil = &(node->n_segs[0]->filaments[0]);
560     else if (dir == E)
561       fil = &(node->e_segs[0]->filaments[0]);
562     else
563       GP_PANIC("add_edge_segs_to_list: bad dir");
564 
565     melem = make_melement(fil->filnumber, fil, signofelem);
566     plist = insert_in_list(melem, plist);
567   }
568 
569   return plist;
570 
571 }
572 
get_or_make_nearest_node(name,index,x,y,z,indsys,gp,node_list)573 NODES *get_or_make_nearest_node(name, index, x, y, z, indsys, gp, node_list)
574      char *name;
575      int index;
576      double x,y,z;
577      SYS *indsys;
578      Nonuni_gp *gp;
579      NPATH *node_list;
580 {
581   G_nodes *nonuni_node;
582   NODES *node;
583   NPATH *nodeL;
584   NODES *makenode();
585 
586   nonuni_node = find_nearest_nonuni_node(x, y, z, gp);
587 
588   if (nonuni_node == NULL) {
589     fprintf(stderr,"Error: Node for plane %s requested at %lg %lg %lg doesn't appear to be within the plane\n",
590 	    gp->grndp->name, x, y, z);
591     exit(1);
592   }
593 
594   /* See if G_nodes node has already been made into a NODES node.
595      If not, creaet it */
596 
597   /* look for nonuni_node already in tempnode list */
598   node = get_nonuni_node_from_list(nonuni_node, node_list);
599 
600   if (node != NULL)
601     /* this non-uni node is the same as another read so far */
602     return node;
603   else {
604     /* make a new one */
605     return make_new_node_with_nonuni(nonuni_node,name, index,
606                                      x, y, z, indsys, gp);
607   }
608 }
609 
610 /* make a new NORMAL node that is derived from a nonuni_node */
make_new_node_with_nonuni(nonuni_node,name,index,x,y,z,indsys,gp)611 NODES *make_new_node_with_nonuni(nonuni_node, name, index, x, y, z, indsys, gp)
612      G_nodes *nonuni_node;
613      char *name;
614      int index;
615      double x,y,z;
616      SYS *indsys;
617      Nonuni_gp *gp;
618 {
619   NODES *node;
620   NODES *makenode();
621 
622   /* make a new one */
623   node = makenode(name, index, x, y, z, NORMAL, indsys);
624   node->gp = gp->grndp;
625   node->gp_node = nonuni_node;
626   return node;
627 }
628 
629 /* find a NODES node in list nodeL whose gp_node is nonuni_node */
630 /* Note: you used to compare getrealnode(nodeL->node)->gp_node != nonuni_node
631    but you took out getrealnode for equiv_rect. Hopefully that won't
632    break anything */
get_nonuni_node_from_list(nonuni_node,nodeL)633 NODES *get_nonuni_node_from_list(nonuni_node, nodeL)
634     G_nodes *nonuni_node;
635     NPATH *nodeL;
636 {
637 
638   while(nodeL != NULL && (nodeL->node)->gp_node != nonuni_node)
639     nodeL = nodeL->next;
640 
641   if (nodeL != NULL)
642     return nodeL->node;
643   else
644     return NULL;
645 }
646 
647 /* return a linked list of all the nodes of the children of
648    this cell that are inside the rectangle with corners (x0,y0),(x1,y1) */
get_nodes_inside_rect(x0,y0,x1,y1,cell,endoflist)649 Llist *get_nodes_inside_rect(x0, y0, x1, y1, cell, endoflist)
650      double x0,y0,x1,y1;
651      Gcell *cell;
652      Llist **endoflist;
653 {
654   Llist *nodelist;
655 
656   /* check if there is any overlap for this child */
657   if (intersection(x0,y0,x1,y1,cell->x0,cell->y0,cell->x1,cell->y1) == FALSE)
658     return NULL;
659 
660   switch(get_children_type(cell)) {
661   case NONE:
662     /* we are a leaf. figure which nodes are inside */
663     return which_nodes_inside(x0,y0,x1,y1,cell, endoflist);
664   case BI:
665     return bi_get_nodes_inside_rect( x0, y0, x1, y1, cell, endoflist);
666     break;
667   case GRID_2D:
668     return grid_get_nodes_inside_rect( x0, y0, x1, y1, cell, endoflist);
669     break;
670   default:
671     GP_PANIC("Unknown child type in set_cell_coords")
672     break;
673   }
674 }
675 
bi_get_nodes_inside_rect(x0,y0,x1,y1,cell,endoflist)676 Llist *bi_get_nodes_inside_rect( x0, y0, x1, y1, cell, endoflist)
677     double x0, y0, x1, y1;
678     Gcell *cell;
679     Llist **endoflist;
680 {
681   Bi *bi = (Bi *)(cell->children);
682   Llist *one;
683   Llist *two;
684   Llist *oneend;
685 
686   one = get_nodes_inside_rect(x0, y0, x1, y1, bi->child1, &oneend);
687   two = get_nodes_inside_rect(x0, y0, x1, y1, bi->child2, endoflist);
688 
689 
690   /* concatenate two lists */
691   if (one == NULL)
692     one = two;
693   else {
694     if (oneend->next != NULL)
695       GP_PANIC("bi_get_nodes_inside_rect: oneend is not end of list!");
696 
697     oneend->next = two;
698 
699     if (two == NULL)
700       *endoflist = oneend;
701   }
702 
703   /* endoflist assigned in two's call */
704   return one;
705 
706 }
707 
grid_get_nodes_inside_rect(x0,y0,x1,y1,cell,endoflist)708 Llist *grid_get_nodes_inside_rect( x0, y0, x1, y1, cell, endoflist)
709     double x0, y0, x1, y1;
710     Gcell *cell;
711     Llist **endoflist;
712 {
713   Grid_2d *grid = (Grid_2d *)cell->children;
714   Llist *current;
715   Llist *current_end;
716   int row_start, row_end, col_start, col_end;
717   int i, j;
718   Llist *entire;
719 
720   *endoflist = NULL;
721 
722   /* to save a little cpu, let's figure out which cells intersect this
723      rectangle */
724 
725   /* get indices that correspond to top left and bottom right corner */
726 
727   get_grid_indices(cell, x0, y1, &row_start, &col_start);  /* top left */
728   get_grid_indices(cell, x1, y0, &row_end, &col_end);  /* bot right */
729 
730   entire = NULL;
731   for(i = row_start; i <= row_end; i++)
732     for(j = col_start; j <= col_end; j++) {
733       current = get_nodes_inside_rect( x0, y0, x1, y1, grid->kids[i][j],
734                                        &current_end);
735       if (current != NULL) {
736         if (*endoflist == NULL)
737           *endoflist = current_end;
738         current_end->next = entire;
739         entire = current;
740       }
741 
742     }
743 
744   return entire;
745 }
746 
747 /* returns the cell indices that contain (x,y) */
get_grid_indices(cell,x,y,pi,pj)748 get_grid_indices(cell, x, y, pi, pj)
749     Gcell *cell;
750     double x,y;
751     int *pi, *pj;
752 {
753   Grid_2d *grid = (Grid_2d *)cell->children;
754   double xn, yn;  /* normalized coords */
755   int y_i;
756 
757   xn = (x - cell->x0)/(cell->x1 - cell->x0);
758   yn = (y - cell->y0)/(cell->y1 - cell->y0);
759 
760   if (xn == 1)
761     *pj = grid->x_cells - 1;
762   else
763     /* round to the nearest cell */
764     *pj = ((int) (xn * grid->x_cells));
765 
766   if (yn == 1)
767     y_i = grid->y_cells - 1;
768   else
769     /* round to the nearest cell */
770     y_i = ((int) (yn * grid->y_cells));
771 
772   *pi = grid->y_cells - y_i - 1;
773 
774 }
775 
776 /* returns TRUE if the two rectangles overlap.
777    lower left=(x0,y0) upper right = (x1,y1) */
intersection(x0,y0,x1,y1,cx0,cy0,cx1,cy1)778 intersection(x0,y0,x1,y1,cx0,cy0,cx1,cy1)
779      double x0,y0,x1,y1,cx0,cy0,cx1,cy1;
780 {
781   if ((x0 < cx0 && x1 < cx0) || (x0 > cx1 && x1 > cx1)
782       && (y0 < cy0 && y1 < cy0) || (y0 > cy1 && y1 > cy1))
783     return FALSE;
784   else
785     return TRUE;
786 
787 }
788 
789 /* returns a linked list of the nodes of this cell contained in the rectangle*/
which_nodes_inside(x0,y0,x1,y1,cell,endoflist)790 Llist *which_nodes_inside(x0,y0,x1,y1,cell, endoflist)
791     double x0, y0, x1, y1;
792     Gcell *cell;
793     Llist **endoflist;
794 {
795   Llist *list = NULL;
796 
797   if (!is_leaf(cell))
798     GP_PANIC("which_nodes_inside: non-leaf cell!");
799 
800   /* is SW inside */
801   if (cell->x0 >= x0 && cell->x0 <= x1 && cell->y0 >= y0 && cell->y0 <= y1) {
802     list = add_ptr_to_list((void *)cell->bndry.nodes[SW], list);
803     if (list->next == NULL)
804       *endoflist = list;
805   }
806 
807   /* is SE inside */
808   if (cell->x1 >= x0 && cell->x1 <= x1 && cell->y0 >= y0 && cell->y0 <= y1) {
809     list = add_ptr_to_list((void *)cell->bndry.nodes[SE], list);
810     if (list->next == NULL)
811       *endoflist = list;
812   }
813 
814   /* is NE inside */
815   if (cell->x1 >= x0 && cell->x1 <= x1 && cell->y1 >= y0 && cell->y1 <= y1) {
816     list = add_ptr_to_list((void *)cell->bndry.nodes[NE], list);
817     if (list->next == NULL)
818       *endoflist = list;
819   }
820 
821   /* is NW inside */
822   if (cell->x0 >= x0 && cell->x0 <= x1 && cell->y1 >= y0 && cell->y1 <= y1) {
823     list = add_ptr_to_list((void *)cell->bndry.nodes[NW], list);
824     if (list->next == NULL)
825       *endoflist = list;
826   }
827 
828   return list;
829 }
830 
free_Llist(list)831 free_Llist(list)
832      Llist *list;
833 {
834   Llist *node;
835 
836   while(list != NULL) {
837     node = list;
838     list = list->next;
839     free(node);
840   }
841 }
842