/* Last edited on 2012-02-25 08:29:21 by stolfilocal */ /* ---------------------------------------------------------------------- */ /* The {clip} option: */ double clip; /* Fraction of paint tray to clip. */ float_image_t *stir_paint(float_image_t *iimg, int NX, int NY, int NW, int ND, double clip, bool_t verbose); float_image_t *oimg = stir_paint(iimg, o->NX, o->NY, o->whorls, o->detailSize, o->clip, o->verbose); /* Map domain of {oimg} to the middle {clip} of the domain of {iimg}: */ double sx = (clip*NXI)/NX; double sy = (clip*NYI)/NY; p->c[0] = 0.5*(1-clip)*NXI + sx*p->c[0]; p->c[1] = 0.5*(1-clip)*NYI + sy*p->c[1]; if (argparser_keyword_present(pp, "-clip")) { o->clip = argparser_get_next_double(pp, 1.0e-6, 1.0); } else { /* Output whole tray: */ o->clip = 1.0; } /* ---------------------------------------------------------------------- */ /* Stirring using {r2_map-expand,r2_map_contract} */ void define_stir_parameters(int NX, int NY, int NW, int ND, r2_t wct[], double wrd[], double wan[], bool_t stretch, bool_t verbose); /* Defines the stirring paramters. The stirring is a series of {NW} whirlpools with decreasing radii and random ceners and angles. Whirlpool number {i} has center {wct[i]}, radius {wrd[i]}, and angle {wan[i]}. The smaller whirlpools have radius {~ND/2}. The center and radius are picked in the domain {[0_NX]×[0_NY]}. If {stretch} is true, the parameters are then mapped with {r2_expand(*,0,NX,0,NY,NULL)}. */ float_image_t *stir_paint(float_image_t *iimg, int NX, int NY, int NW, int ND, double clip, bool_t verbose) { bool_t debug = TRUE; bool_t stretch = FALSE; int NC = iimg->sz[0]; int NXI = iimg->sz[1]; int NYI = iimg->sz[2]; float_image_t *oimg = float_image_new(NC, NX, NY); ix_reduction_t red = (stretch ? ix_reduction_SINGLE : ix_reduction_MIRROR); bool_t avg = TRUE; int order = 0; /* Interpolation order */ /* Stirring parameters: */ r2_t wct[NW]; /* Whirlpool centers. */ double wrd[NW]; /* Whirlpool centers. */ double wan[NW]; /* Whirlpool angles. */ define_stir_parameters(NXI, NYI, NW, ND, wct, wrd, wan, stretch, verbose); auto void stir_map (r2_t *p, r2x2_t *J); /* Stores in {q} the point of {[0_NX]×[0_NY]} that ends up at {p} after the stirring .*/ float_image_transform_all(iimg, red, &stir_map, 0.0, avg, order, NULL, oimg); if (debug) { bool_t round = TRUE; bool_t diagonal = FALSE; int sub = 4; int drad = 2.0; /* Dot radius. */ int dhwd = 0.5; /* Pen half-width. */ int i; for (i = 0; i < NW; i++) { int c; r2_t ctr = wct[i]; double rad = wrd[i]; if (stretch) { /* Convert {ctr,rad} from stretched to image coords: */ r2_t dir = ctr; r2_dir(&dir, &dir); r2_t ref; r2_mix(1.0, &ctr, rad, &dir, &ref); r2_map_contract(&(ctr), 0, NX, 0, NY, NULL); r2_map_contract(&(ref), 0, NX, 0, NY, NULL); rad = r2_dist(&ctr, &ref); } for (c = 0; c < NC; c++) { float_image_paint_dot(oimg, c, ctr.c[0], ctr.c[1], drad, 0, round, diagonal, 1.0, NAN, sub); float_image_paint_dot(oimg, c, ctr.c[0], ctr.c[1], rad, dhwd, round, diagonal, NAN, 1.0, sub); } } } return oimg; void stir_map (r2_t *p, r2x2_t *J) { /* Map domain of {oimg} to the middle {clip} of the domain of {iimg}: */ double sx = (clip*NXI)/NX; double sy = (clip*NYI)/NY; p->c[0] = 0.5*(1-clip)*NXI + sx*p->c[0]; p->c[1] = 0.5*(1-clip)*NYI + sy*p->c[1]; if (J != NULL) { r2x2_t PJ; r2x2_ident(&PJ); PJ.c[0][0] = sx; PJ.c[1][1] = sy; r2x2_mul(J, &PJ, J); } if (stretch) { r2_map_expand(p, 0, NXI, 0, NYI, J); } int i; for (i = 0; i < NW; i++) { r2_map_twirl(p, &(wct[i]), wrd[i], wan[i], J); } if (stretch) { r2_map_contract(p, 0, NXI, 0, NYI, J); } } } void define_stir_parameters(int NX, int NY, int NW, int ND, r2_t wct[], double wrd[], double wan[], bool_t stretch, bool_t verbose) { r2_t ctrImg = (r2_t){{ 0.5*NX, 0.5*NY }}; /* Image center. */ double radMax = 0.20*fmin(NX, NY); /* Max whorl radius. */ double radIni = 0.20*hypot(NX, NY); /* Ideal radius of first whorl. */ double radFin = 0.50*fmin(ND,radIni); /* Ideal radius of last whorl. */ double scaFin = radFin/radIni; /* Ideal scale from last to first. */ assert(scaFin < 1.0); double offset = (scaFin*(NW-1))/(1 - scaFin); /* Offset for denominator. */ int i, k; for (i = 0; i < NW; i++) { /* Decide the relative scale of this whorl: */ double scale = offset/(i + offset); /* Pick a turn angle {ang} in radians: */ double ang = 0; for (k = 0; k < 4; k++) { double ank = sqrt(scale)*3*M_PI*(drandom()-drandom()); if (fabs(ank) > fabs(ang)) { ang = ank; } } /* Pick a nominal radius {rad} in image space: */ double rad = fmin(radMax, scale*radIni); r2_t ctr; if (stretch) { double mrg = fmin(rad, 0.40*fmin(NX,NY)); assert(2*mrg < NX); assert(2*mrg < NY); /* Pick the center {ctr} in image space: */ double ctrx = 0.5*(drandom()+drandom()); double ctry = 0.5*(drandom()+drandom()); ctr = (r2_t){{ mrg + (NX-2*mrg)*ctrx, mrg + (NY-2*mrg)*ctry }}; if (verbose) { r2_gen_print(stderr, &ctr, "%14.8f", "raw ctr = ( ", " ", " )"); fprintf(stderr, " rad = %14.8f ang = %+12.8f (%+8.2f°)\n", rad, ang, ang*180/M_PI); } /* Get a point {ref} on the {ctr,rad}-circle, towards {ctrImg}, in image space: */ r2_t dir; r2_sub(&ctrImg, &ctr, &dir); r2_dir(&dir, &dir); r2_t ref; r2_mix(1.0, &ctr, rad, &dir, &ref); /* Now map {ctr,ref} to expanded space ({R^2}): */ r2_map_expand(&(ctr), 0, NX, 0, NY, NULL); r2_map_expand(&(ref), 0, NX, 0, NY, NULL); /* Get the radius {rad} in expanded space: */ rad = r2_dist(&ctr, &ref); if (verbose) { r2_gen_print(stderr, &ctr, "%+14.6f", "fin ctr = ( ", " ", " )"); fprintf(stderr, " rad = %14.8f ang = %+12.8f (%+8.2f°)\n", rad, ang, ang*180/M_PI); } } else { ctr = (r2_t){{ NX*drandom(),NY*drandom() }}; if (verbose) { r2_gen_print(stderr, &ctr, "%14.8f", "fin ctr = ( ", " ", " )"); fprintf(stderr, " rad = %14.8f ang = %+12.8f (%+8.2f°)\n", rad, ang, ang*180/M_PI); } } /* Set param vector: */ wct[i] = ctr; wrd[i] = rad; wan[i] = ang; } }