| 9 | 
  | 
 */ | 
| 10 | 
  | 
 | 
| 11 | 
  | 
#include "glare.h" | 
| 12 | 
+ | 
#include "linregr.h" | 
| 13 | 
  | 
 | 
| 14 | 
  | 
#define vcont(vd,vu)    ((vu)-(vd)<=SEPS) | 
| 15 | 
  | 
#define hcont(s1,s2)    ((s1)->r-(s2)->l>=-SEPS&&(s2)->r-(s1)->l>=-SEPS) | 
| 51 | 
  | 
        for (v = vsize; v >= -vsize; v--) { | 
| 52 | 
  | 
                close_sources(v); | 
| 53 | 
  | 
#ifndef DEBUG | 
| 54 | 
< | 
                if (verbose) | 
| 54 | 
> | 
                if (verbose) { | 
| 55 | 
  | 
                        fprintf(stderr, "%s: analyzing... %3ld%%\r", | 
| 56 | 
  | 
                                progname, 100L*(vsize-v)/(2*vsize)); | 
| 57 | 
+ | 
                        fflush(stderr); | 
| 58 | 
+ | 
                } | 
| 59 | 
  | 
#endif | 
| 60 | 
  | 
                getviewspan(v, spanbr); | 
| 61 | 
  | 
                left = hsize + 1; | 
| 278 | 
  | 
} | 
| 279 | 
  | 
 | 
| 280 | 
  | 
 | 
| 281 | 
+ | 
struct srcspan * | 
| 282 | 
+ | 
splitspan(sso, h, v, m)         /* divide source span at point */ | 
| 283 | 
+ | 
register struct srcspan *sso; | 
| 284 | 
+ | 
double  h, v, m; | 
| 285 | 
+ | 
{ | 
| 286 | 
+ | 
        register struct srcspan *ssn; | 
| 287 | 
+ | 
        double  d; | 
| 288 | 
+ | 
        int     hs; | 
| 289 | 
+ | 
 | 
| 290 | 
+ | 
        d = h - m*(sso->v - v); | 
| 291 | 
+ | 
        hs = d < 0. ? d-.5 : d+.5; | 
| 292 | 
+ | 
        if (sso->l >= hs) | 
| 293 | 
+ | 
                return(NULL); | 
| 294 | 
+ | 
        if (sso->r <= hs) | 
| 295 | 
+ | 
                return(sso); | 
| 296 | 
+ | 
                                /* need to split it */ | 
| 297 | 
+ | 
        ssn = (struct srcspan *)malloc(sizeof(struct srcspan)); | 
| 298 | 
+ | 
        if (ssn == NULL) | 
| 299 | 
+ | 
                memerr("source spans in splitspan"); | 
| 300 | 
+ | 
        ssn->brsum = (double)(hs - sso->l)/(sso->r - sso->l) * sso->brsum; | 
| 301 | 
+ | 
        sso->brsum -= ssn->brsum; | 
| 302 | 
+ | 
        ssn->v = sso->v; | 
| 303 | 
+ | 
        ssn->l = sso->l; | 
| 304 | 
+ | 
        ssn->r = sso->l = hs; | 
| 305 | 
+ | 
        return(ssn); | 
| 306 | 
+ | 
} | 
| 307 | 
+ | 
 | 
| 308 | 
+ | 
 | 
| 309 | 
+ | 
struct source * | 
| 310 | 
+ | 
splitsource(so)                 /* divide source in two if it's big and long */ | 
| 311 | 
+ | 
struct source   *so; | 
| 312 | 
+ | 
{ | 
| 313 | 
+ | 
        LRSUM   lr; | 
| 314 | 
+ | 
        LRLIN   fit; | 
| 315 | 
+ | 
        register struct srcspan *ss, *ssn; | 
| 316 | 
+ | 
        struct srcspan  *ssl, *ssnl, head; | 
| 317 | 
+ | 
        int     h; | 
| 318 | 
+ | 
        double  mh, mv; | 
| 319 | 
+ | 
        struct source   *sn; | 
| 320 | 
+ | 
 | 
| 321 | 
+ | 
        lrclear(&lr); | 
| 322 | 
+ | 
        for (ss = so->first; ss != NULL; ss = ss->next) | 
| 323 | 
+ | 
                for (h = ss->l; h < ss->r; h++) | 
| 324 | 
+ | 
                        lrpoint(h, ss->v, &lr); | 
| 325 | 
+ | 
        if ((double)lr.n/(sampdens*sampdens) < SABIG) | 
| 326 | 
+ | 
                return(NULL);                   /* too small */ | 
| 327 | 
+ | 
        if (lrfit(&fit, &lr) < 0) | 
| 328 | 
+ | 
                return(NULL);                   /* can't fit a line */ | 
| 329 | 
+ | 
        if (fit.correlation < LCORR && fit.correlation > -LCORR) | 
| 330 | 
+ | 
                return(NULL); | 
| 331 | 
+ | 
        if (verbose) | 
| 332 | 
+ | 
                fprintf(stderr, "%s: splitting large source\n", progname); | 
| 333 | 
+ | 
        mh = lrxavg(&lr); | 
| 334 | 
+ | 
        mv = lryavg(&lr); | 
| 335 | 
+ | 
        sn = (struct source *)malloc(sizeof(struct source)); | 
| 336 | 
+ | 
        if (sn == NULL) | 
| 337 | 
+ | 
                memerr("source records in splitsource"); | 
| 338 | 
+ | 
        sn->dom = 0.0; | 
| 339 | 
+ | 
        sn->first = NULL; | 
| 340 | 
+ | 
        ssnl = NULL; | 
| 341 | 
+ | 
        head.next = so->first; | 
| 342 | 
+ | 
        ssl = &head; | 
| 343 | 
+ | 
        for (ss = so->first; ss != NULL; ssl = ss, ss = ss->next) | 
| 344 | 
+ | 
                if ((ssn = splitspan(ss, mh, mv, fit.slope)) != NULL) { | 
| 345 | 
+ | 
                        if (ssn == ss) {        /* remove from old */ | 
| 346 | 
+ | 
                                ssl->next = ss->next; | 
| 347 | 
+ | 
                                ss = ssl; | 
| 348 | 
+ | 
                        } | 
| 349 | 
+ | 
                        if (ssnl == NULL)       /* add to new */ | 
| 350 | 
+ | 
                                sn->first = ssn; | 
| 351 | 
+ | 
                        else | 
| 352 | 
+ | 
                                ssnl->next = ssn; | 
| 353 | 
+ | 
                        ssn->next = NULL; | 
| 354 | 
+ | 
                        ssnl = ssn; | 
| 355 | 
+ | 
                } | 
| 356 | 
+ | 
        so->first = head.next; | 
| 357 | 
+ | 
        return(sn); | 
| 358 | 
+ | 
} | 
| 359 | 
+ | 
 | 
| 360 | 
+ | 
 | 
| 361 | 
  | 
donesource(sp)                  /* finished with this source */ | 
| 362 | 
  | 
register struct source  *sp; | 
| 363 | 
  | 
{ | 
| 364 | 
+ | 
        struct source   *newsrc; | 
| 365 | 
  | 
        register struct srcspan *ss; | 
| 366 | 
  | 
        int     h, n; | 
| 367 | 
  | 
        double  hsum, vsum, d; | 
| 368 | 
  | 
 | 
| 369 | 
+ | 
        while ((newsrc = splitsource(sp)) != NULL)      /* split it? */ | 
| 370 | 
+ | 
                donesource(newsrc); | 
| 371 | 
  | 
        sp->dom = 0.0; | 
| 372 | 
  | 
        hsum = vsum = 0.0; | 
| 373 | 
  | 
        sp->brt = 0.0; | 
| 381 | 
  | 
                        vsum += d*ss->v; | 
| 382 | 
  | 
                        sp->dom += d; | 
| 383 | 
  | 
                } | 
| 298 | 
– | 
                free((char *)ss); | 
| 384 | 
  | 
        } | 
| 385 | 
< | 
        sp->first = NULL; | 
| 385 | 
> | 
        freespans(sp); | 
| 386 | 
> | 
        if (sp->dom <= FTINY) {         /* must be right at edge of image */ | 
| 387 | 
> | 
                free((char *)sp); | 
| 388 | 
> | 
                return; | 
| 389 | 
> | 
        } | 
| 390 | 
  | 
        sp->brt /= (double)n; | 
| 391 | 
  | 
        compdir(sp->dir, (int)(hsum/sp->dom), (int)(vsum/sp->dom)); | 
| 392 | 
  | 
        sp->next = donelist; | 
| 457 | 
  | 
                indirect[i].sum += d * s->brt; | 
| 458 | 
  | 
                indirect[i].n += d; | 
| 459 | 
  | 
        } | 
| 460 | 
< | 
        for ( ; s->first != NULL; s->first = s->first->next) | 
| 372 | 
< | 
                free((char *)s->first); | 
| 460 | 
> | 
        freespans(s); | 
| 461 | 
  | 
        free((char *)s); | 
| 462 | 
+ | 
} | 
| 463 | 
+ | 
 | 
| 464 | 
+ | 
 | 
| 465 | 
+ | 
freespans(sp)                   /* free spans associated with source */ | 
| 466 | 
+ | 
struct source   *sp; | 
| 467 | 
+ | 
{ | 
| 468 | 
+ | 
        register struct srcspan *ss; | 
| 469 | 
+ | 
 | 
| 470 | 
+ | 
        while ((ss = sp->first) != NULL) { | 
| 471 | 
+ | 
                sp->first = ss->next; | 
| 472 | 
+ | 
                free((char *)ss); | 
| 473 | 
+ | 
        } | 
| 474 | 
  | 
} |