Piece-wise radial uniformly random distribution of points around a reference point

Question

I have been learning Asymptote, which is a vector graphics programming language. In addition to its nice integration with LaTeX and friends, it is capable of generating 3D rendered artworks. For a semi-learning semi-research purpose, I'd like to generate some random points around a center point, called preMigrationCenter such that their distribution vary in different radial ranges from that point. In particular, for two radial ranges [r1,r2] and [r2,r3], if r1 < r2 < r3, then the density of points spread in [r1,r2] shall be higher than that of [r2,r3], with respect to some hyperparameter settings. The code below is the kernel of a very larger file which is the reason that some hyperparemeters, e.g., the overall size of the page, may seem larger than what it should be in this simplified case. Except the general stylistic comments on my code, I shall be grateful if one can share any optimization ideas on the random generation process of my code, say, the way I manage seed points, if possible.

General specifications:

TeX engine and unit size are set.

settings.tex="pdflatex";
unitsize(1mm);

Moreover, using the code borrowed from here, absolute page size is set to alleviate any future scaling issues.

pair page_min = (0, 0);
pair page_max = (500, 250);
clip(box(page_min, page_max));
fixedscaling(page_min, page_max);

Hyperparameters:

Here are the parameters to set a population for my points, radii, and the proportion of the sub-population associated with each radial range. Points base and preMigrationCenter are the reference point of the whole schematic and the center to define radii, respectively.

real totalRadius = 150;
real mostInnerRadius = totalRadius*0.1;
real innerRadius = totalRadius*0.2;
real outerRadius = totalRadius*0.3;
real mostOuterRadius = totalRadius*0.4;

real mostInnerPopulation = 30;
real innerPopulation = 25;
real outerPopulation = 20;
real mostOuterPopulation = 15;

pair base = (0,0);
pair preMigrationCenter = (150, 125);

Utility functions:

Each point has three properties, collectively known as configuration, all of which have to be set in a random manner, including its coordinate components and its orientation. The function below, given a radial range [lowerBound, upperBound], uses randunit() function of Asymptote. unitrand() returns a random number sampled from a uniform distribution whose range is [0,1]. The output of this function is the radius component of the polar coordinate of the desired random point.

real getRandomRadius(real lowerBound, real upperBound, int i)
{
    srand(seconds()+i);
    real randomRadius = (upperBound - lowerBound)*unitrand() + lowerBound; 
    return randomRadius;
}

The angle component of the desired point is computed by

real getRandomTheta(int i)
{
    srand(seconds()+11*(i+3));
    real randomTheta = 360*unitrand();
    return randomTheta;
}

The code below randomly generates an orientation for the point, transforms its coordinate to a Cartesian frame, and returns a packaged configuration. (A point obviously has no orientation but, as one sees later, each point hosts a png image which is not symmetric, thereby requiring a specified orientation.)

real[] computeConfiguration(pair migrationCenter, real randomRadius, real randomTheta, int i)
{
    srand(seconds()+i*19);
    real xCoord = migrationCenter.x + randomRadius*Cos(randomTheta);
    real yCoord = migrationCenter.y + randomRadius*Sin(randomTheta);
    real orientation = unitrand()*360;
    real[] configuration = {xCoord, yCoord, orientation};
    return configuration;
}

Here is a plotter. The image shape is first scaled, and then rotated according to the orientation of its configuration.

void plot(real[] configuration)
{   
    label(rotate(configuration[2], base)*graphic("shape.png","scale=0.03"), (configuration[0], configuration[1]));
}

Finally, the wrapper below manages the whole process for a given radial range.

void process(pair migrationCenter, real marginalPopulation, real smallerRadius, real largerRadius)
{
    for (int i = 1; i < marginalPopulation + 1; ++i)
    {
        real randomRadius = getRandomRadius(smallerRadius, largerRadius, i);
        real randomTheta = getRandomTheta(i);
        real[] configuration = computeConfiguration(migrationCenter, randomRadius, randomTheta, i);
        plot(configuration);
    }
}

Operations:

In the end, we call as many process({arguments}) functions as the number of radial ranges. sleep() functions are embedded to assure that the seed points generated by srand functions are different.

path preMostOuterCircle = circle(preMigrationCenter, mostOuterRadius);
draw(preMostOuterCircle, red+dashed);

process(preMigrationCenter, mostInnerPopulation, 0, mostInnerRadius);
sleep(1);
process(preMigrationCenter, innerPopulation, mostInnerRadius, innerRadius);
sleep(2);
process(preMigrationCenter, outerPopulation, innerRadius, outerRadius);
sleep(1);
process(preMigrationCenter, mostOuterPopulation, outerRadius, mostOuterRadius);

Execution:

Gluing all snippets above together in <file_name>.asy,

///////////////////////////  General specifications /////////////////////////

settings.tex="pdflatex";
unitsize(1mm);

pair page_min = (0, 0);
pair page_max = (500, 250);
clip(box(page_min, page_max));
fixedscaling(page_min, page_max);

///////////////////////////// Hyperparameters /////////////////////////////

real totalRadius = 150;
real mostInnerRadius = totalRadius*0.1;
real innerRadius = totalRadius*0.2;
real outerRadius = totalRadius*0.3;
real mostOuterRadius = totalRadius*0.4;

real mostInnerPopulation = 30;
real innerPopulation = 25;
real outerPopulation = 20;
real mostOuterPopulation = 15;

pair base = (0,0);
pair preMigrationCenter = (150, 125);

//////////////////////////// Utility functions /////////////////////////// 

real getRandomRadius(real lowerBound, real upperBound, int i)
{
    srand(seconds()+i);
    real randomRadius = (upperBound - lowerBound)*unitrand() + lowerBound; 
    return randomRadius;
}

real getRandomTheta(int i)
{
    srand(seconds()+11*(i+3));
    real randomTheta = 360*unitrand();
    return randomTheta;
}

real[] computeConfiguration(pair migrationCenter, real randomRadius, real randomTheta, int i)
{
    srand(seconds()+i*19);
    real xCoord = migrationCenter.x + randomRadius*Cos(randomTheta);
    real yCoord = migrationCenter.y + randomRadius*Sin(randomTheta);
    real orientation = unitrand()*360;
    real[] configuration = {xCoord, yCoord, orientation};
    return configuration;
}

void plot(real[] configuration)
{   
    label(rotate(configuration[2], base)*graphic("ant.png","scale=0.03"), (configuration[0], configuration[1]));
}

void process(pair migrationCenter, real marginalPopulation, real smallerRadius, real largerRadius)
{
    for (int i = 1; i < marginalPopulation + 1; ++i)
    {
        real randomRadius = getRandomRadius(smallerRadius, largerRadius, i);
        real randomTheta = getRandomTheta(i);
        real[] configuration = computeConfiguration(migrationCenter, randomRadius, randomTheta, i);
        plot(configuration);
    }
}

/////////////////////////////// Operations /////////////////////////////// 

path preMostOuterCircle = circle(preMigrationCenter, mostOuterRadius);
draw(preMostOuterCircle, red+dashed);

process(preMigrationCenter, mostInnerPopulation, 0, mostInnerRadius);
sleep(1);
process(preMigrationCenter, innerPopulation, mostInnerRadius, innerRadius);
sleep(2);
process(preMigrationCenter, outerPopulation, innerRadius, outerRadius);
sleep(1);
process(preMigrationCenter, mostOuterPopulation, outerRadius, mostOuterRadius);

one may execute

asy.exe -f pdf <file_name>.asy

to get something like

Any comments on the code quality, readibility, and performance, especially the way I have managed the random generations using multiple srand() calls, is appreciated.

Answer 1

    srand(seconds()+i);

    srand(seconds()+11*(i+3));

    srand(seconds()+i*19);

The way that pseudorandom generators normally work is that you seed them (sometimes implicitly) and then they produce a predictable but changing sequence.

Here though, if you call these functions with the same value of i and have not crossed a second boundary, you will get the same "random" value. Presumably that's why you added the extra stuff like +11*(i+3) so as to at least get three different values.

A better approach (although still not cryptographically secure) is to call srand once. I would put it at the beginning of the section that you call "Operations". A simple

srand(seconds());
path preMostOuterCircle = circle(preMigrationCenter, mostOuterRadius);

The second line is unchanged from your original code and only included for context.

Then you can get rid of the rest of the srand calls and the sleep calls.

A more complicated solution would define some kind of flag variable globally as false. Then each of these functions could check if the flag is true. If not, they could call srand(seconds()); and set the flag to true. It might perhaps be possible to do that in an Asymptote package that would include the random functions (getRandomRadius, getRandomTheta, and computeConfiguration).

It might also be possible to call srand automatically when a package is imported. That would be both simple and robust if you only call rand/unitrand from one package.

Whether you do it by manually invoking srand as part of your operations or automatically in a package, what's important is that you call it once per program invocation before ever calling unitrand (or rand). Note that you can define the functions first. What's important is when it is invoked for the first time. It doesn't matter when the function is defined; it does matter when it is called.

If you want to avoid being able to invoke the code twice quickly with the same results, you might consider instead using srand(cputime().change.clock));. But still call it once per program invocation.

I'm not seeing documentation of how Asymptote implements srand. So I'm guessing that it basically calls the C or C++ srand. It seems to be implemented in C++ (GitHub)