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 ¤t_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