@param node the node to begin computing from @param i @param z @param workSpace the indicies are the accumulated contribution to the gradient sans multiplicative terms in the first 2 indices. @return the contribution to the normalizing constant Z
(Quadtree.Node node, int i, double[] z, double[] workSpace)
| 525 | * @return the contribution to the normalizing constant Z |
| 526 | */ |
| 527 | private double computeF_rep(Quadtree.Node node, int i, double[] z, double[] workSpace) |
| 528 | { |
| 529 | if(node == null || node.N_cell == 0 || node.indx == i) |
| 530 | return 0; |
| 531 | /* |
| 532 | * Original paper says to use the diagonal divided by the squared 2 |
| 533 | * norm. This dosn't seem to work at all. Tried some different ideas |
| 534 | * with 0.5 as the threshold until I found one that worked. |
| 535 | * Squaring the values would normally not be helpful, but since we are working with tiny values it makes them smaller, making it easier to hit the go |
| 536 | */ |
| 537 | double x = z[i*2]; |
| 538 | double y = z[i*2+1]; |
| 539 | // double r_cell = node.diagLen(); |
| 540 | double r_cell = Math.max(node.maxX-node.minX, node.maxY-node.minY); |
| 541 | r_cell*=r_cell; |
| 542 | double mass_x = node.x_mass/node.N_cell; |
| 543 | double mass_y = node.y_mass/node.N_cell; |
| 544 | double dot = (mass_x-x)*(mass_x-x)+(mass_y-y)*(mass_y-y); |
| 545 | |
| 546 | |
| 547 | if(node.NW == null || r_cell < theta*dot)//good enough! |
| 548 | { |
| 549 | if(node.indx == i) |
| 550 | return 0; |
| 551 | |
| 552 | double Z = 1.0/(1.0 + dot); |
| 553 | double q_cell_Z_sqrd = -node.N_cell*(Z*Z); |
| 554 | |
| 555 | workSpace[0] += q_cell_Z_sqrd*(x-mass_x); |
| 556 | workSpace[1] += q_cell_Z_sqrd*(y-mass_y); |
| 557 | return Z*node.N_cell; |
| 558 | } |
| 559 | else//further subdivide |
| 560 | { |
| 561 | double Z_sum = 0; |
| 562 | for(Quadtree.Node child : node) |
| 563 | Z_sum += computeF_rep(child, i, z, workSpace); |
| 564 | return Z_sum; |
| 565 | } |
| 566 | } |
| 567 | |
| 568 | /** |
| 569 | * |