slices.pd (PDL-2.074) | : | slices.pd (PDL-2.075) | ||
---|---|---|---|---|
skipping to change at line 98 | skipping to change at line 98 | |||
C<index>, C<index1d>, and C<index2d> provide rudimentary index indirection. | C<index>, C<index1d>, and C<index2d> provide rudimentary index indirection. | |||
=for example | =for example | |||
$c = index($source,$ind); | $c = index($source,$ind); | |||
$c = index1d($source,$ind); | $c = index1d($source,$ind); | |||
$c = index2d($source2,$ind1,$ind2); | $c = index2d($source2,$ind1,$ind2); | |||
use the C<$ind> variables as indices to look up values in C<$source>. | use the C<$ind> variables as indices to look up values in C<$source>. | |||
The three routines thread slightly differently. | The three routines broadcast slightly differently. | |||
=over 3 | =over 3 | |||
=item * | =item * | |||
C<index> uses direct threading for 1-D indexing across the 0 dim | C<index> uses direct broadcasting for 1-D indexing across the 0 dim | |||
of C<$source>. It can thread over source thread dims or index thread | of C<$source>. It can broadcast over source broadcast dims or index broadcast | |||
dims, but not (easily) both: If C<$source> has more than 1 | dims, but not (easily) both: If C<$source> has more than 1 | |||
dimension and C<$ind> has more than 0 dimensions, they must agree in | dimension and C<$ind> has more than 0 dimensions, they must agree in | |||
a threading sense. | a broadcasting sense. | |||
=item * | =item * | |||
C<index1d> uses a single active dim in C<$ind> to produce a list of | C<index1d> uses a single active dim in C<$ind> to produce a list of | |||
indexed values in the 0 dim of the output - it is useful for | indexed values in the 0 dim of the output - it is useful for | |||
collapsing C<$source> by indexing with a single row of values along | collapsing C<$source> by indexing with a single row of values along | |||
C<$source>'s 0 dimension. The output has the same number of dims as | C<$source>'s 0 dimension. The output has the same number of dims as | |||
C<$source>. The 0 dim of the output has size 1 if C<$ind> is a | C<$source>. The 0 dim of the output has size 1 if C<$ind> is a | |||
scalar, and the same size as the 0 dim of C<$ind> if it is not. If | scalar, and the same size as the 0 dim of C<$ind> if it is not. If | |||
C<$ind> and C<$source> both have more than 1 dim, then all dims higher | C<$ind> and C<$source> both have more than 1 dim, then all dims higher | |||
than 0 must agree in a threading sense. | than 0 must agree in a broadcasting sense. | |||
=item * | =item * | |||
C<index2d> works like C<index> but uses separate ndarrays for X and Y | C<index2d> works like C<index> but uses separate ndarrays for X and Y | |||
coordinates. For more general N-dimensional indexing, see the | coordinates. For more general N-dimensional indexing, see the | |||
L<PDL::NiceSlice> syntax or L<PDL::Slices> (in particular C<slice>, | L<PDL::NiceSlice> syntax or L<PDL::Slices> (in particular C<slice>, | |||
C<indexND>, and C<range>). | C<indexND>, and C<range>). | |||
=back | =back | |||
These functions are two-way, i.e. after | These functions are two-way, i.e. after | |||
$c = $x->index(pdl[0,5,8]); | $c = $x->index(pdl[0,5,8]); | |||
$c .= pdl [0,2,4]; | $c .= pdl [0,2,4]; | |||
the changes in C<$c> will flow back to C<$x>. | the changes in C<$c> will flow back to C<$x>. | |||
C<index> provids simple threading: multiple-dimensioned arrays are treated | C<index> provids simple broadcasting: multiple-dimensioned arrays are treated | |||
as collections of 1-D arrays, so that | as collections of 1-D arrays, so that | |||
$x = xvals(10,10)+10*yvals(10,10); | $x = xvals(10,10)+10*yvals(10,10); | |||
$y = $x->index(3); | $y = $x->index(3); | |||
$c = $x->index(9-xvals(10)); | $c = $x->index(9-xvals(10)); | |||
puts a single column from C<$x> into C<$y>, and puts a single element | puts a single column from C<$x> into C<$y>, and puts a single element | |||
from each column of C<$x> into C<$c>. If you want to extract multiple | from each column of C<$x> into C<$c>. If you want to extract multiple | |||
columns from an array in one operation, see L</dice> or | columns from an array in one operation, see L</dice> or | |||
L</indexND>. | L</indexND>. | |||
skipping to change at line 510 | skipping to change at line 510 | |||
example, the four expressions [0,1], ['forbid','truncate'], ['f','t'], | example, the four expressions [0,1], ['forbid','truncate'], ['f','t'], | |||
and 'ft' all specify that violating the boundary in the 0th dimension | and 'ft' all specify that violating the boundary in the 0th dimension | |||
throws an error, and all other dimensions get truncated. | throws an error, and all other dimensions get truncated. | |||
If you feed in a single string, it is interpreted as a packed boundary | If you feed in a single string, it is interpreted as a packed boundary | |||
array if all of its characters are valid boundary specifiers (e.g. 'pet'), | array if all of its characters are valid boundary specifiers (e.g. 'pet'), | |||
but as a single word-style specifier if they are not (e.g. 'forbid'). | but as a single word-style specifier if they are not (e.g. 'forbid'). | |||
B<OUTPUT> | B<OUTPUT> | |||
The output threads over both C<$index> and C<$source>. Because implicit | The output broadcasts over both C<$index> and C<$source>. Because implicit | |||
threading can happen in a couple of ways, a little thought is needed. The | broadcasting can happen in a couple of ways, a little thought is needed. The | |||
returned dimension list is stacked up like this: | returned dimension list is stacked up like this: | |||
(index thread dims), (index dims (size)), (source thread dims) | (index broadcast dims), (index dims (size)), (source broadcast dims) | |||
The first few dims of the output correspond to the extra dims of | The first few dims of the output correspond to the extra dims of | |||
C<$index> (beyond the 0 dim). They allow you to pick out individual | C<$index> (beyond the 0 dim). They allow you to pick out individual | |||
ranges from a large, threaded collection. | ranges from a large, broadcasted collection. | |||
The middle few dims of the output correspond to the size dims | The middle few dims of the output correspond to the size dims | |||
specified in C<$size>, and contain the range of values that is extracted | specified in C<$size>, and contain the range of values that is extracted | |||
at each location in C<$source>. Every nonzero element of C<$size> is copied to | at each location in C<$source>. Every nonzero element of C<$size> is copied to | |||
the dimension list here, so that if you feed in (for example) C<$size | the dimension list here, so that if you feed in (for example) C<$size | |||
= [2,0,1]> you get an index dim list of C<(2,1)>. | = [2,0,1]> you get an index dim list of C<(2,1)>. | |||
The last few dims of the output correspond to extra dims of C<$source> beyond | The last few dims of the output correspond to extra dims of C<$source> beyond | |||
the number of dims indexed by C<$index>. These dims act like ordinary | the number of dims indexed by C<$index>. These dims act like ordinary | |||
thread dims, because adding more dims to C<$source> just tacks extra dims | broadcast dims, because adding more dims to C<$source> just tacks extra dims | |||
on the end of the output. Each source thread dim ranges over the entire | on the end of the output. Each source broadcast dim ranges over the entire | |||
corresponding dim of C<$source>. | corresponding dim of C<$source>. | |||
B<Dataflow>: Dataflow is bidirectional. | B<Dataflow>: Dataflow is bidirectional. | |||
B<Examples>: | B<Examples>: | |||
Here are basic examples of C<range> operation, showing how to get | Here are basic examples of C<range> operation, showing how to get | |||
ranges out of a small matrix. The first few examples show extraction | ranges out of a small matrix. The first few examples show extraction | |||
and selection of individual chunks. The last example shows | and selection of individual chunks. The last example shows | |||
how to mark loci in the original matrix (using dataflow). | how to mark loci in the original matrix (using dataflow). | |||
skipping to change at line 580 | skipping to change at line 580 | |||
[11 22] | [11 22] | |||
[23 1] | [23 1] | |||
] | ] | |||
[ | [ | |||
[21 32] | [21 32] | |||
[33 11] | [33 11] | |||
] | ] | |||
] | ] | |||
] | ] | |||
pdl> $src = xvals(5,3)*10+yvals(5,3) | pdl> $src = xvals(5,3)*10+yvals(5,3) | |||
pdl> print $src->range(3,1) # Thread over y dimension in $src | pdl> print $src->range(3,1) # Broadcast over y dimension in $src | |||
[ | [ | |||
[30] | [30] | |||
[31] | [31] | |||
[32] | [32] | |||
] | ] | |||
pdl> $src = zeroes(5,4); | pdl> $src = zeroes(5,4); | |||
pdl> $src->range(pdl([2,3],[0,1]),pdl(2,1)) .= xvals(2,2,1) + 1 | pdl> $src->range(pdl([2,3],[0,1]),pdl(2,1)) .= xvals(2,2,1) + 1 | |||
pdl> print $src | pdl> print $src | |||
[ | [ | |||
skipping to change at line 654 | skipping to change at line 654 | |||
$x = $data->range($index, $sizes, $bound)->sever; | $x = $data->range($index, $sizes, $bound)->sever; | |||
$aa = $data->range($index, $sizes, $bound); | $aa = $data->range($index, $sizes, $bound); | |||
map { $x($_ - 1) .= $_; } (1..$x->nelem); # Lots of little ops | map { $x($_ - 1) .= $_; } (1..$x->nelem); # Lots of little ops | |||
$aa .= $x; | $aa .= $x; | |||
C<range> is a perl front-end to a PP function, C<rangeb>. Calling | C<range> is a perl front-end to a PP function, C<rangeb>. Calling | |||
C<rangeb> is marginally faster but requires that you include all arguments. | C<rangeb> is marginally faster but requires that you include all arguments. | |||
DEVEL NOTES | DEVEL NOTES | |||
* index thread dimensions are effectively clumped internally. This | * index broadcast dimensions are effectively clumped internally. This | |||
makes it easier to loop over the index array but a little more brain-bending | makes it easier to loop over the index array but a little more brain-bending | |||
to tease out the algorithm. | to tease out the algorithm. | |||
=cut | =cut | |||
EOD | EOD | |||
HandleBad => 1, | HandleBad => 1, | |||
TwoWay => 1, | TwoWay => 1, | |||
P2Child => 1, | P2Child => 1, | |||
# | # | |||
# rdim: dimensionality of each range (0 dim of index PDL) | # rdim: dimensionality of each range (0 dim of index PDL) | |||
# | # | |||
# ntsize: number of nonzero size dimensions | # ntsize: number of nonzero size dimensions | |||
# sizes: array of range sizes, indexed (0..rdim-1). A zero element means | # sizes: array of range sizes, indexed (0..rdim-1). A zero element means | |||
# that the dimension is omitted from the child dim list. | # that the dimension is omitted from the child dim list. | |||
# corners: parent coordinates of each corner, running fastest over coord index. | # corners: parent coordinates of each corner, running fastest over coord index. | |||
# (indexed 0 .. (nitems-1)*(rdim)+rdim-1) | # (indexed 0 .. (nitems-1)*(rdim)+rdim-1) | |||
# nitems: total number of list elements (product of itdims) | # nitems: total number of list elements (product of itdims) | |||
# itdim: number of index thread dimensions | # itdim: number of index broadcast dimensions | |||
# itdims: Size of each index thread dimension, indexed (0..itdim-1) | # itdims: Size of each index broadcast dimension, indexed (0..itdim-1) | |||
# | # | |||
# bsize: Number of independently specified boundary conditions | # bsize: Number of independently specified boundary conditions | |||
# nsizes: Number of independently specified range dim sizes | # nsizes: Number of independently specified range dim sizes | |||
# boundary: Array containing all the boundary condition specs | # boundary: Array containing all the boundary condition specs | |||
# indord: Order/size of the indexing dim (0th dim of $index) | # indord: Order/size of the indexing dim (0th dim of $index) | |||
Comp => 'PDL_Indx rdim; | Comp => 'PDL_Indx rdim; | |||
PDL_Indx nitems; | PDL_Indx nitems; | |||
PDL_Indx itdim; | PDL_Indx itdim; | |||
PDL_Indx ntsize; | PDL_Indx ntsize; | |||
PDL_Indx bsize; | PDL_Indx bsize; | |||
PDL_Indx nsizes; | PDL_Indx nsizes; | |||
PDL_Indx sizes[$COMP(rdim)]; | PDL_Indx sizes[$COMP(rdim)]; | |||
PDL_Indx itdims[$COMP(itdim)]; | PDL_Indx itdims[$COMP(itdim)]; | |||
PDL_Indx corners[$COMP(rdim) * $COMP(nitems)]; | PDL_Indx corners[$COMP(rdim) * $COMP(nitems)]; | |||
char boundary[$COMP(rdim)]; | char boundary[$COMP(rdim)]; | |||
char size_pdl_destroy; | ||||
char ind_pdl_destroy; | ||||
', | ', | |||
MakeComp => <<'EOD-MakeComp', | MakeComp => <<'EOD-MakeComp', | |||
pdl *size_pdl; | pdl *size_pdl; | |||
PDL_RETERROR(PDL_err, PDL->make_physdims(ind_pdl)); | PDL_RETERROR(PDL_err, PDL->make_physdims(ind_pdl)); | |||
$COMP(size_pdl_destroy) = $COMP(ind_pdl_destroy) = 0; | ||||
/* Generalized empties are ok, but not in the special 0 dim (the index vector) * / | /* Generalized empties are ok, but not in the special 0 dim (the index vector) * / | |||
if(ind_pdl->dims[0] == 0) | if(ind_pdl->dims[0] == 0) | |||
{ $CROAK("can't handle Empty indices -- call range instead"); } | { $CROAK("can't handle Empty indices -- call range instead"); } | |||
/*** | /*** | |||
* Ensure that the index is a PDL_Indx. If there's no loss of information, | * Ensure that the index is a PDL_Indx. If there's no loss of information, | |||
* just upgrade it -- otherwise, make a temporary copy. | * just upgrade it -- otherwise, make a temporary copy. | |||
*/ | */ | |||
switch(ind_pdl->datatype) { | switch(ind_pdl->datatype) { | |||
default: /* Most types: */ | default: /* Most types: */ | |||
ind_pdl = PDL->hard_copy(ind_pdl); /* copy and fall through */ | ind_pdl = PDL->hard_copy(ind_pdl); /* copy and fall through */ | |||
if (!ind_pdl) $CROAK("Error in hard_copy"); | if (!ind_pdl) $CROAK("Error in hard_copy"); | |||
$COMP(ind_pdl_destroy) = 1; | ||||
case PDL_SB: case PDL_B: case PDL_S: case PDL_US: case PDL_L: case PDL_UL: case PDL_LL: case PDL_ULL: | case PDL_SB: case PDL_B: case PDL_S: case PDL_US: case PDL_L: case PDL_UL: case PDL_LL: case PDL_ULL: | |||
PDL_RETERROR(PDL_err, PDL->converttype(ind_pdl,PDL_IND)); /* convert in place . */ | PDL_RETERROR(PDL_err, PDL->converttype(ind_pdl,PDL_IND)); /* convert in place . */ | |||
break; | break; | |||
case PDL_IND: | case PDL_IND: | |||
PDL_RETERROR(PDL_err, PDL->make_physical(ind_pdl)); | PDL_RETERROR(PDL_err, PDL->make_physical(ind_pdl)); | |||
break; | break; | |||
} | } | |||
/*** | /*** | |||
* Figure sizes of the COMP arrrays and allocate them. | * Figure sizes of the COMP arrrays and allocate them. | |||
skipping to change at line 776 | skipping to change at line 780 | |||
* first time it gets hit. | * first time it gets hit. | |||
*/ | */ | |||
$CROAK("Unknown boundary condition '%c' in range",bstr[i]); | $CROAK("Unknown boundary condition '%c' in range",bstr[i]); | |||
} | } | |||
break; | break; | |||
} // end of switch | } // end of switch | |||
} | } | |||
} | } | |||
} | } | |||
/*** | /*** | |||
* Store the sizes of the index-thread dims | * Store the sizes of the index-broadcast dims | |||
*/ | */ | |||
{ | { | |||
PDL_Indx i; | PDL_Indx i; | |||
PDL_Indx nd = ind_pdl->ndims - 1; | PDL_Indx nd = ind_pdl->ndims - 1; | |||
for(i=0; i < nd ; i++) | for(i=0; i < nd ; i++) | |||
$COMP(itdims[i]) = ind_pdl->dims[i+1]; | $COMP(itdims[i]) = ind_pdl->dims[i+1]; | |||
} | } | |||
/*** | /*** | |||
* Check and condition the size ndarray, and store sizes of the ranges | * Check and condition the size ndarray, and store sizes of the ranges | |||
skipping to change at line 815 | skipping to change at line 819 | |||
for(i=0;i<$COMP(rdim);i++) | for(i=0;i<$COMP(rdim);i++) | |||
$COMP(sizes[i]) = 0; | $COMP(sizes[i]) = 0; | |||
} else { | } else { | |||
// Convert size PDL to PDL_IND to support indices | // Convert size PDL to PDL_IND to support indices | |||
switch(size_pdl->datatype) { | switch(size_pdl->datatype) { | |||
default: /* Most types: */ | default: /* Most types: */ | |||
size_pdl = PDL->hard_copy(size_pdl); /* copy and fall through */ | size_pdl = PDL->hard_copy(size_pdl); /* copy and fall through */ | |||
if (!size_pdl) $CROAK("Error in hard_copy"); | if (!size_pdl) $CROAK("Error in hard_copy"); | |||
$COMP(size_pdl_destroy) = 1; | ||||
case PDL_SB: case PDL_B: case PDL_S: case PDL_US: case PDL_L: case PDL_UL: case PDL_LL: case PDL_ULL: | case PDL_SB: case PDL_B: case PDL_S: case PDL_US: case PDL_L: case PDL_UL: case PDL_LL: case PDL_ULL: | |||
PDL_RETERROR(PDL_err, PDL->converttype(size_pdl,PDL_IND)); /* convert in place. */ | PDL_RETERROR(PDL_err, PDL->converttype(size_pdl,PDL_IND)); /* convert in place. */ | |||
break; | break; | |||
case PDL_IND: | case PDL_IND: | |||
PDL_RETERROR(PDL_err, PDL->make_physical(size_pdl)); | PDL_RETERROR(PDL_err, PDL->make_physical(size_pdl)); | |||
break; | break; | |||
} | } | |||
$COMP(nsizes) = size_pdl->nvals; /* Store for later permissiveness check * / | $COMP(nsizes) = size_pdl->nvals; /* Store for later permissiveness check * / | |||
skipping to change at line 866 | skipping to change at line 871 | |||
if($COMP(sizes[i])) | if($COMP(sizes[i])) | |||
ntsize++; | ntsize++; | |||
$COMP(ntsize) = ntsize; | $COMP(ntsize) = ntsize; | |||
} | } | |||
/*** | /*** | |||
* Stash coordinates of the corners | * Stash coordinates of the corners | |||
*/ | */ | |||
PDL_Indx j,k,ioff, iter[$COMP(itdim)]; | PDL_Indx j,k,ioff, iter[$COMP(itdim)]; | |||
/* initialize iterator to loop over index threads */ | /* initialize iterator to loop over index broadcasts */ | |||
PDL_Indx *cptr = iter; | PDL_Indx *cptr = iter; | |||
for(k=0;k<$COMP(itdim);k++) | for(k=0;k<$COMP(itdim);k++) | |||
*(cptr++) = 0; | *(cptr++) = 0; | |||
cptr = $COMP(corners); | cptr = $COMP(corners); | |||
do { | do { | |||
/* accumulate offset into the index from the iterator */ | /* accumulate offset into the index from the iterator */ | |||
for(k=ioff=0;k<$COMP(itdim);k++) | for(k=ioff=0;k<$COMP(itdim);k++) | |||
ioff += iter[k] * ind_pdl->dimincs[k+1]; | ioff += iter[k] * ind_pdl->dimincs[k+1]; | |||
/* Loop over the 0th dim of index, copying coords. */ | /* Loop over the 0th dim of index, copying coords. */ | |||
/* This is the natural place to check for permissive ranging; too */ | /* This is the natural place to check for permissive ranging; too */ | |||
/* bad we don't have access to the parent ndarray here... */ | /* bad we don't have access to the parent ndarray here... */ | |||
for(j=0;j<$COMP(rdim);j++) | for(j=0;j<$COMP(rdim);j++) | |||
*(cptr++) = ((PDL_Indx *)(ind_pdl->data))[ioff + ind_pdl->dimincs[0] * j]; | *(cptr++) = ((PDL_Indx *)(ind_pdl->data))[ioff + ind_pdl->dimincs[0] * j]; | |||
/* Increment the iterator -- the test increments, the body carries. */ | /* Increment the iterator -- the test increments, the body carries. */ | |||
for(k=0; k<$COMP(itdim) && (++(iter[k]))>=($COMP(itdims)[k]) ;k++) | for(k=0; k<$COMP(itdim) && (++(iter[k]))>=($COMP(itdims)[k]) ;k++) | |||
iter[k] = 0; | iter[k] = 0; | |||
} while(k<$COMP(itdim)); | } while(k<$COMP(itdim)); | |||
if ($COMP(ind_pdl_destroy)) | ||||
PDL->destroy(ind_pdl); /* finished with our copy */ | ||||
if ($COMP(size_pdl_destroy)) | ||||
PDL->destroy(size_pdl); /* finished with our copy */ | ||||
EOD-MakeComp | EOD-MakeComp | |||
RedoDims => <<'EOD-RedoDims' , | RedoDims => <<'EOD-RedoDims' , | |||
PDL_Indx stdim = $PARENT(ndims) - $COMP(rdim); | PDL_Indx stdim = $PARENT(ndims) - $COMP(rdim); | |||
PDL_Indx dim,inc; | PDL_Indx dim,inc; | |||
PDL_Indx i,rdvalid; | PDL_Indx i,rdvalid; | |||
// Speed bump for ludicrous cases | // Speed bump for ludicrous cases | |||
if( $COMP(rdim) > $PARENT(ndims)+5 && $COMP(nsizes) != $COMP(rdim)) { | if( $COMP(rdim) > $PARENT(ndims)+5 && $COMP(nsizes) != $COMP(rdim)) { | |||
$CROAK("Ludicrous number of extra dims in range index; leaving child null.\n (%d implicit dims is > 5; index has %d dims; source has %d dim%s.)\n This often means that your index PDL is incorrect. To avoid this message,\n allo cate dummy dims in the source or use %d dims in range's size field.\n",$COMP(rdi m)-$PARENT(ndims),$COMP(rdim),$PARENT(ndims),($PARENT(ndims))>1?"s":"",$COMP(rdi m)); | $CROAK("Ludicrous number of extra dims in range index; leaving child null.\n (%d implicit dims is > 5; index has %d dims; source has %d dim%s.)\n This often means that your index PDL is incorrect. To avoid this message,\n allo cate dummy dims in the source or use %d dims in range's size field.\n",$COMP(rdi m)-$PARENT(ndims),$COMP(rdim),$PARENT(ndims),($PARENT(ndims))>1?"s":"",$COMP(rdi m)); | |||
skipping to change at line 914 | skipping to change at line 923 | |||
/* Copy size dimensions to child, crunching as we go. */ | /* Copy size dimensions to child, crunching as we go. */ | |||
dim = $COMP(itdim); | dim = $COMP(itdim); | |||
for(i=rdvalid=0;i<$COMP(rdim);i++) { | for(i=rdvalid=0;i<$COMP(rdim);i++) { | |||
if($COMP(sizes[i])) { | if($COMP(sizes[i])) { | |||
rdvalid++; | rdvalid++; | |||
$CHILD(dimincs[dim]) = inc; | $CHILD(dimincs[dim]) = inc; | |||
inc *= ($CHILD(dims[dim++]) = $COMP(sizes[i])); /* assignment */ | inc *= ($CHILD(dims[dim++]) = $COMP(sizes[i])); /* assignment */ | |||
} | } | |||
} | } | |||
/* Copy index thread dimensions to child */ | /* Copy index broadcast dimensions to child */ | |||
for(dim=0; dim<$COMP(itdim); dim++) { | for(dim=0; dim<$COMP(itdim); dim++) { | |||
$CHILD(dimincs[dim]) = inc; | $CHILD(dimincs[dim]) = inc; | |||
inc *= ($CHILD(dims[dim]) = $COMP(itdims[dim])); /* assignment */ | inc *= ($CHILD(dims[dim]) = $COMP(itdims[dim])); /* assignment */ | |||
} | } | |||
/* Copy source thread dimensions to child */ | /* Copy source broadcast dimensions to child */ | |||
dim = $COMP(itdim) + rdvalid; | dim = $COMP(itdim) + rdvalid; | |||
for(i=0;i<stdim;i++) { | for(i=0;i<stdim;i++) { | |||
$CHILD(dimincs[dim]) = inc; | $CHILD(dimincs[dim]) = inc; | |||
inc *= ($CHILD(dims[dim++]) = $PARENT(dims[i+$COMP(rdim)])); /* assignment * / | inc *= ($CHILD(dims[dim++]) = $PARENT(dims[i+$COMP(rdim)])); /* assignment * / | |||
} | } | |||
/* Cover bizarre case where the source PDL is empty - in that case, change */ | /* Cover bizarre case where the source PDL is empty - in that case, change */ | |||
/* all non-barfing boundary conditions to truncation, since we have no data */ | /* all non-barfing boundary conditions to truncation, since we have no data */ | |||
/* to reflect, extend, or mirror. */ | /* to reflect, extend, or mirror. */ | |||
if($PARENT(dims[0]) == 0) { | if($PARENT(dims[0]) == 0) { | |||
skipping to change at line 949 | skipping to change at line 958 | |||
$SETDIMS(); | $SETDIMS(); | |||
EOD-RedoDims | EOD-RedoDims | |||
EquivCPOffsCode => <<'EOD-EquivCPOffsCode', | EquivCPOffsCode => <<'EOD-EquivCPOffsCode', | |||
PDL_Indx *ip; /* vector iterator */ | PDL_Indx *ip; /* vector iterator */ | |||
PDL_Indx *sp; /* size vector including stdims */ | PDL_Indx *sp; /* size vector including stdims */ | |||
PDL_Indx *coords; /* current coordinates */ | PDL_Indx *coords; /* current coordinates */ | |||
PDL_Indx k; /* index */ | PDL_Indx k; /* index */ | |||
PDL_Indx item; /* index thread iterator */ | PDL_Indx item; /* index broadcast iterator */ | |||
PDL_Indx pdim = $PDL(PARENT)->ndims; | PDL_Indx pdim = $PDL(PARENT)->ndims; | |||
PDL_Indx rdim = $COMP(rdim); | PDL_Indx rdim = $COMP(rdim); | |||
PDL_Indx prdim = PDLMIN(rdim,pdim); | PDL_Indx prdim = PDLMIN(rdim,pdim); | |||
PDL_Indx iter2[pdim * 2 + rdim]; | PDL_Indx iter2[pdim * 2 + rdim]; | |||
PDL_Indx *sizes = iter2 + pdim; | PDL_Indx *sizes = iter2 + pdim; | |||
coords = sizes + pdim; | coords = sizes + pdim; | |||
/* Figure out size vector */ | /* Figure out size vector */ | |||
for(ip = $COMP(sizes), sp = sizes, k=0; k<rdim; k++) | for(ip = $COMP(sizes), sp = sizes, k=0; k<rdim; k++) | |||
*(sp++) = *(ip++); | *(sp++) = *(ip++); | |||
skipping to change at line 1057 | skipping to change at line 1066 | |||
/* dims caused by permissive ranging. */ | /* dims caused by permissive ranging. */ | |||
coff = $PDL(CHILD)->dimincs[0] * item; | coff = $PDL(CHILD)->dimincs[0] * item; | |||
for(k2 = $COMP(itdim), poff = k = 0; | for(k2 = $COMP(itdim), poff = k = 0; | |||
k < prdim; | k < prdim; | |||
k++) { | k++) { | |||
poff += coords[k]*$PDL(PARENT)->dimincs[k]; | poff += coords[k]*$PDL(PARENT)->dimincs[k]; | |||
if($COMP(sizes[k])) | if($COMP(sizes[k])) | |||
coff += iter2[k] * $PDL(CHILD)->dimincs[k2++]; | coff += iter2[k] * $PDL(CHILD)->dimincs[k2++]; | |||
} | } | |||
/* Loop the copy over all the source thread dims (above rdim). */ | /* Loop the copy over all the source broadcast dims (above rdim). */ | |||
do { | do { | |||
PDL_Indx poff1 = poff; | PDL_Indx poff1 = poff; | |||
PDL_Indx coff1 = coff; | PDL_Indx coff1 = coff; | |||
/* Accumulate the offset due to source threading */ | /* Accumulate the offset due to source broadcasting */ | |||
for(k2 = $COMP(itdim) + $COMP(ntsize), k = rdim; | for(k2 = $COMP(itdim) + $COMP(ntsize), k = rdim; | |||
k < pdim; | k < pdim; | |||
k++) { | k++) { | |||
poff1 += iter2[k] * $PDL(PARENT)->dimincs[k]; | poff1 += iter2[k] * $PDL(PARENT)->dimincs[k]; | |||
coff1 += iter2[k] * $PDL(CHILD)->dimincs[k2++]; | coff1 += iter2[k] * $PDL(CHILD)->dimincs[k2++]; | |||
} | } | |||
/* Finally -- make the copy | /* Finally -- make the copy | |||
* EQUIVCPTRUNC works like EQUIVCPOFFS but with checking for | * EQUIVCPTRUNC works like EQUIVCPOFFS but with checking for | |||
* out-of-bounds conditions. | * out-of-bounds conditions. | |||
*/ | */ | |||
$EQUIVCPTRUNC(coff1,poff1,trunc); | $EQUIVCPTRUNC(coff1,poff1,trunc); | |||
/* Increment the source thread iterator */ | /* Increment the source broadcast iterator */ | |||
for( k=$COMP(rdim); | for( k=$COMP(rdim); | |||
k < pdim && (++(iter2[k]) >= $PDL(PARENT)->dims[k]); | k < pdim && (++(iter2[k]) >= $PDL(PARENT)->dims[k]); | |||
k++) | k++) | |||
iter2[k] = 0; | iter2[k] = 0; | |||
} while(k < pdim); /* end of source-thread iteration */ | } while(k < pdim); /* end of source-broadcast iteration */ | |||
/* Increment the in-range iterator */ | /* Increment the in-range iterator */ | |||
for(k = 0; | for(k = 0; | |||
k < $COMP(rdim) && (++(iter2[k]) >= $COMP(sizes[k])); | k < $COMP(rdim) && (++(iter2[k]) >= $COMP(sizes[k])); | |||
k++) | k++) | |||
iter2[k] = 0; | iter2[k] = 0; | |||
} while(k < $COMP(rdim)); /* end of main iteration */ | } while(k < $COMP(rdim)); /* end of main iteration */ | |||
} /* end of item do loop */ | } /* end of item do loop */ | |||
EOD-EquivCPOffsCode | EOD-EquivCPOffsCode | |||
skipping to change at line 1109 | skipping to change at line 1118 | |||
'rld', | 'rld', | |||
GenericTypes => [ppdefs_all], | GenericTypes => [ppdefs_all], | |||
Pars=>'indx a(n); b(n); [o]c(m);', | Pars=>'indx a(n); b(n); [o]c(m);', | |||
PMCode =><<'EOD', | PMCode =><<'EOD', | |||
sub PDL::rld { | sub PDL::rld { | |||
my ($x,$y) = @_; | my ($x,$y) = @_; | |||
my ($c); | my ($c); | |||
if ($#_ == 2) { | if ($#_ == 2) { | |||
$c = $_[2]; | $c = $_[2]; | |||
} else { | } else { | |||
# XXX Need to improve emulation of threading in auto-generating c | # XXX Need to improve emulation of broadcasting in auto-generating c | |||
my ($size) = $x->sumover->max->sclr; | my ($size) = $x->sumover->max->sclr; | |||
my (@dims) = $x->dims; | my (@dims) = $x->dims; | |||
shift @dims; | shift @dims; | |||
$c = $y->zeroes($size,@dims); | $c = $y->zeroes($size,@dims); | |||
} | } | |||
&PDL::_rld_int($x,$y,$c); | &PDL::_rld_int($x,$y,$c); | |||
$c; | $c; | |||
} | } | |||
EOD | EOD | |||
Code=>' | Code=>' | |||
skipping to change at line 1154 | skipping to change at line 1163 | |||
EOD | EOD | |||
); | ); | |||
=head2 rle | =head2 rle | |||
=cut | =cut | |||
pp_def( | pp_def( | |||
'rle', | 'rle', | |||
GenericTypes => [ppdefs_all], | GenericTypes => [ppdefs_all], | |||
Pars=>'c(n); indx [o]a(m); [o]b(m);', | Pars=>'c(n); indx [o]a(m); [o]b(m);', | |||
#this RedoDimsCode sets $SIZE(m)==$SIZE(n), but the slice in the PMCode below ma | RedoDimsCode=>'$SIZE(m)=$SIZE(n);', | |||
kes m<=n. | ||||
RedoDimsCode=>'$SIZE(m)=$PDL(c)->dims[0];', | ||||
PMCode=><<'EOC', | PMCode=><<'EOC', | |||
sub PDL::rle { | sub PDL::rle { | |||
my $c = shift; | my $c = shift; | |||
my ($x,$y) = @_==2 ? @_ : (null,null); | my ($x,$y) = @_==2 ? @_ : (null,null); | |||
PDL::_rle_int($c,$x,$y); | PDL::_rle_int($c,$x,$y); | |||
my $max_ind = ($c->ndims<2) ? ($x!=0)->sumover-1 : | my $max_ind = ($c->ndims<2) ? ($x!=0)->sumover-1 : | |||
($x!=0)->clump(1..$x->ndims-1)->sumover->max->sc lr-1; | ($x!=0)->clump(1..$x->ndims-1)->sumover->max->sc lr-1; | |||
return ($x->slice("0:$max_ind"),$y->slice("0:$max_ind")); | return ($x->slice("0:$max_ind"),$y->slice("0:$max_ind")); | |||
} | } | |||
EOC | EOC | |||
skipping to change at line 1196 | skipping to change at line 1204 | |||
', | ', | |||
Doc => <<'EOD' | Doc => <<'EOD' | |||
=for ref | =for ref | |||
Run-length encode a vector | Run-length encode a vector | |||
Given vector C<$c>, generate a vector C<$x> with the number of each | Given vector C<$c>, generate a vector C<$x> with the number of each | |||
element, and a vector C<$y> of the unique values. New in PDL 2.017, | element, and a vector C<$y> of the unique values. New in PDL 2.017, | |||
only the elements up to the first instance of C<0> in C<$x> are | only the elements up to the first instance of C<0> in C<$x> are | |||
returned, which makes the common use case of a 1-dimensional C<$c> simpler. | returned, which makes the common use case of a 1-dimensional C<$c> simpler. | |||
For threaded operation, C<$x> and C<$y> will be large enough | For broadcast operation, C<$x> and C<$y> will be large enough | |||
to hold the largest row of C<$y>, and only the elements up to the | to hold the largest row of C<$y>, and only the elements up to the | |||
first instance of C<0> in each row of C<$x> should be considered. | first instance of C<0> in each row of C<$x> should be considered. | |||
=for example | =for example | |||
$c = floor(4*random(10)); | $c = floor(4*random(10)); | |||
rle($c,$x=null,$y=null); | rle($c,$x=null,$y=null); | |||
#or | #or | |||
($x,$y) = rle($c); | ($x,$y) = rle($c); | |||
skipping to change at line 1245 | skipping to change at line 1253 | |||
P2Child => 1, | P2Child => 1, | |||
RedoDims => 'int i; PDL_Indx d1; | RedoDims => 'int i; PDL_Indx d1; | |||
/* truncate overly long clumps to just clump existing dimensions */ | /* truncate overly long clumps to just clump existing dimensions */ | |||
if($COMP(n) > $PARENT(ndims)) | if($COMP(n) > $PARENT(ndims)) | |||
$COMP(n) = $PARENT(ndims); | $COMP(n) = $PARENT(ndims); | |||
if($COMP(n) < -1) | if($COMP(n) < -1) | |||
$COMP(n) = $PARENT(ndims) + $COMP(n) + 1; | $COMP(n) = $PARENT(ndims) + $COMP(n) + 1; | |||
PDL_Indx nrem = ($COMP(n) == -1 ? $PARENT(threadids[0]) : $COMP (n)); | PDL_Indx nrem = ($COMP(n) == -1 ? $PARENT(broadcastids[0]) : $C OMP(n)); | |||
$SETNDIMS($PARENT(ndims) - nrem + 1); | $SETNDIMS($PARENT(ndims) - nrem + 1); | |||
d1=1; | d1=1; | |||
for(i=0; i<nrem; i++) { | for(i=0; i<nrem; i++) { | |||
d1 *= $PARENT(dims[i]); | d1 *= $PARENT(dims[i]); | |||
} | } | |||
$CHILD(dims[0]) = d1; | $CHILD(dims[0]) = d1; | |||
for(; i<$PARENT(ndims); i++) { | for(; i<$PARENT(ndims); i++) { | |||
$CHILD(dims[i-nrem+1]) = $PARENT(dims[i]); | $CHILD(dims[i-nrem+1]) = $PARENT(dims[i]); | |||
} | } | |||
$SETDIMS(); | $SETDIMS(); | |||
$SETDELTATHREADIDS(1-nrem); | $SETDELTABROADCASTIDS(1-nrem); | |||
', | ', | |||
EquivCPOffsCode => ' | EquivCPOffsCode => ' | |||
PDL_Indx i; | PDL_Indx i; | |||
for(i=0; i<$PDL(CHILD)->nvals; i++) { | for(i=0; i<$PDL(CHILD)->nvals; i++) { | |||
$EQUIVCPOFFS(i,i); | $EQUIVCPOFFS(i,i); | |||
} | } | |||
', | ', | |||
TwoWay => 1, | TwoWay => 1, | |||
Doc => 'internal', | Doc => 'internal', | |||
); | ); | |||
skipping to change at line 1278 | skipping to change at line 1286 | |||
=head2 xchg | =head2 xchg | |||
=cut | =cut | |||
pp_def( | pp_def( | |||
'xchg', | 'xchg', | |||
OtherPars => 'PDL_Indx n1; PDL_Indx n2;', | OtherPars => 'PDL_Indx n1; PDL_Indx n2;', | |||
TwoWay => 1, | TwoWay => 1, | |||
P2Child => 1, | P2Child => 1, | |||
AffinePriv => 1, | AffinePriv => 1, | |||
EquivDimCheck => 'if ($COMP(n1) <0) | EquivDimCheck => 'if ($COMP(n1) <0) | |||
$COMP(n1) += $PARENT(threadids[0]); | $COMP(n1) += $PARENT(broadcastids[0]); | |||
if ($COMP(n2) <0) | if ($COMP(n2) <0) | |||
$COMP(n2) += $PARENT(threadids[0]); | $COMP(n2) += $PARENT(broadcastids[0]); | |||
if (PDLMIN($COMP(n1),$COMP(n2)) <0 || | if (PDLMIN($COMP(n1),$COMP(n2)) <0 || | |||
PDLMAX($COMP(n1),$COMP(n2)) >= $PARENT(threadids[0])) | PDLMAX($COMP(n1),$COMP(n2)) >= $PARENT(broadcastids[0])) | |||
$CROAK("One of dims %d, %d out of range: should be 0<=dim<%d", | $CROAK("One of dims %d, %d out of range: should be 0<=dim<%d", | |||
$COMP(n1),$COMP(n2),$PARENT(threadids[0]));', | $COMP(n1),$COMP(n2),$PARENT(broadcastids[0]));', | |||
EquivPDimExpr => '(($CDIM == $COMP(n1)) ? $COMP(n2) : ($CDIM == $COMP(n2 )) ? $COMP(n1) : $CDIM)', | EquivPDimExpr => '(($CDIM == $COMP(n1)) ? $COMP(n2) : ($CDIM == $COMP(n2 )) ? $COMP(n1) : $CDIM)', | |||
Doc => <<'EOD', | Doc => <<'EOD', | |||
=for ref | =for ref | |||
exchange two dimensions | exchange two dimensions | |||
Negative dimension indices count from the end. | Negative dimension indices count from the end. | |||
The command | The command | |||
skipping to change at line 1424 | skipping to change at line 1432 | |||
# Check to make sure all the dims are within bounds | # Check to make sure all the dims are within bounds | |||
for my $i(0..$#newDimOrder) { | for my $i(0..$#newDimOrder) { | |||
my $dim = $newDimOrder[$i]; | my $dim = $newDimOrder[$i]; | |||
if($dim < 0 || $dim > $#newDimOrder) { | if($dim < 0 || $dim > $#newDimOrder) { | |||
my $errString = "PDL::reorder: Dim index $newDimOrder[$i] out of r ange in position $i\n(range is 0-$#newDimOrder)"; | my $errString = "PDL::reorder: Dim index $newDimOrder[$i] out of r ange in position $i\n(range is 0-$#newDimOrder)"; | |||
barf($errString); | barf($errString); | |||
} | } | |||
} | } | |||
# Checking that they are all present and also not duplicated is done by thread() [I think] | # Checking that they are all present and also not duplicated is done by broadcast() [I think] | |||
# a quicker way to do the reorder | # a quicker way to do the reorder | |||
return $pdl->thread(@newDimOrder)->unthread(0); | return $pdl->broadcast(@newDimOrder)->unbroadcast(0); | |||
} | } | |||
EOD | EOD | |||
=head2 mv | =head2 mv | |||
=cut | =cut | |||
pp_addhdr(<<'EOF'); | pp_addhdr(<<'EOF'); | |||
#define EQUIVDIM(dima,dimb,cdim,inc) \ | #define EQUIVDIM(dima,dimb,cdim,inc) \ | |||
((cdim < PDLMIN(dima,dimb) || cdim > PDLMAX(dima,dimb)) ? \ | ((cdim < PDLMIN(dima,dimb) || cdim > PDLMAX(dima,dimb)) ? \ | |||
cdim : ((cdim == dimb) ? dima : cdim + inc)) | cdim : ((cdim == dimb) ? dima : cdim + inc)) | |||
EOF | EOF | |||
pp_def( | pp_def( | |||
'mv', | 'mv', | |||
OtherPars => 'PDL_Indx n1; PDL_Indx n2;', | OtherPars => 'PDL_Indx n1; PDL_Indx n2;', | |||
TwoWay => 1, | TwoWay => 1, | |||
P2Child => 1, | P2Child => 1, | |||
AffinePriv => 1, | AffinePriv => 1, | |||
EquivDimCheck => 'if ($COMP(n1) <0) | EquivDimCheck => 'if ($COMP(n1) <0) | |||
$COMP(n1) += $PARENT(threadids[0]); | $COMP(n1) += $PARENT(broadcastids[0]); | |||
if ($COMP(n2) <0) | if ($COMP(n2) <0) | |||
$COMP(n2) += $PARENT(threadids[0]); | $COMP(n2) += $PARENT(broadcastids[0]); | |||
if (PDLMIN($COMP(n1),$COMP(n2)) <0 || | if (PDLMIN($COMP(n1),$COMP(n2)) <0 || | |||
PDLMAX($COMP(n1),$COMP(n2)) >= $PARENT(threadids[0])) | PDLMAX($COMP(n1),$COMP(n2)) >= $PARENT(broadcastids[0])) | |||
$CROAK("One of dims %d, %d out of range: should be 0<=dim<%d", | $CROAK("One of dims %d, %d out of range: should be 0<=dim<%d", | |||
$COMP(n1),$COMP(n2),$PARENT(threadids[0]));', | $COMP(n1),$COMP(n2),$PARENT(broadcastids[0]));', | |||
EquivPDimExpr => '( | EquivPDimExpr => '( | |||
$COMP(n1) == $COMP(n2) ? $CDIM : | $COMP(n1) == $COMP(n2) ? $CDIM : | |||
$COMP(n1) < $COMP(n2) ? EQUIVDIM($COMP(n1),$COMP(n2),$CDIM,1) : | $COMP(n1) < $COMP(n2) ? EQUIVDIM($COMP(n1),$COMP(n2),$CDIM,1) : | |||
EQUIVDIM($COMP(n1),$COMP(n2),$CDIM,-1) | EQUIVDIM($COMP(n1),$COMP(n2),$CDIM,-1) | |||
)', | )', | |||
Doc => << 'EOD', | Doc => << 'EOD', | |||
=for ref | =for ref | |||
move a dimension to another position | move a dimension to another position | |||
skipping to change at line 1603 | skipping to change at line 1611 | |||
$d = $x->diagonal(dim1, dim2,...) | $d = $x->diagonal(dim1, dim2,...) | |||
=for example | =for example | |||
$y = $x->diagonal(0,2,5); | $y = $x->diagonal(0,2,5); | |||
the ndarray C<$y> has dimensions C<(5,3,4,6)> and | the ndarray C<$y> has dimensions C<(5,3,4,6)> and | |||
C<$y-E<gt>at(2,1,0,1)> refers | C<$y-E<gt>at(2,1,0,1)> refers | |||
to C<$x-E<gt>at(2,1,2,0,1,2)>. | to C<$x-E<gt>at(2,1,2,0,1,2)>. | |||
NOTE: diagonal doesn't handle threadids correctly. XXX FIX | NOTE: diagonal doesn't handle broadcastids correctly. XXX FIX | |||
pdl> $x = zeroes(3,3,3); | pdl> $x = zeroes(3,3,3); | |||
pdl> ($y = $x->diagonal(0,1))++; | pdl> ($y = $x->diagonal(0,1))++; | |||
pdl> p $x | pdl> p $x | |||
[ | [ | |||
[ | [ | |||
[1 0 0] | [1 0 0] | |||
[0 1 0] | [0 1 0] | |||
[0 0 1] | [0 0 1] | |||
] | ] | |||
skipping to change at line 1798 | skipping to change at line 1806 | |||
$y(n=>j) = $x(n=>i); | $y(n=>j) = $x(n=>i); | |||
}', | }', | |||
BackCode=>$rotate_code.' | BackCode=>$rotate_code.' | |||
$x(n=>i) = $y(n=>j); | $x(n=>i) = $y(n=>j); | |||
} | } | |||
' | ' | |||
); | ); | |||
# This is a bit tricky. Hope I haven't missed any cases. | # This is a bit tricky. Hope I haven't missed any cases. | |||
pp_def( | pp_def( | |||
'threadI', | 'broadcastI', | |||
Doc => <<'EOD', | Doc => <<'EOD', | |||
=for ref | =for ref | |||
internal | internal | |||
Put some dimensions to a threadid. | Put some dimensions to a broadcastid. | |||
=for example | =for example | |||
$y = $x->threadI(0,1,5); # thread over dims 1,5 in id 1 | $y = $x->broadcastI(0,1,5); # broadcast over dims 1,5 in id 1 | |||
=cut | =cut | |||
EOD | EOD | |||
P2Child => 1, | P2Child => 1, | |||
TwoWay => 1, | TwoWay => 1, | |||
AffinePriv => 1, | AffinePriv => 1, | |||
OtherPars => "PDL_Indx id; PDL_Indx whichdims[]", | OtherPars => "PDL_Indx id; PDL_Indx whichdims[]", | |||
Comp => 'PDL_Indx nrealwhichdims', | Comp => 'PDL_Indx nrealwhichdims', | |||
MakeComp => ' | MakeComp => ' | |||
skipping to change at line 1841 | skipping to change at line 1849 | |||
} | } | |||
', | ', | |||
RedoDims => ' | RedoDims => ' | |||
PDL_Indx nthc,i,j,flag; | PDL_Indx nthc,i,j,flag; | |||
$SETNDIMS($PARENT(ndims)); | $SETNDIMS($PARENT(ndims)); | |||
$DOPRIVALLOC(); | $DOPRIVALLOC(); | |||
$PRIV(offs) = 0; | $PRIV(offs) = 0; | |||
nthc=0; | nthc=0; | |||
for(i=0; i<$PARENT(ndims); i++) { | for(i=0; i<$PARENT(ndims); i++) { | |||
flag=0; | flag=0; | |||
if($PARENT(nthreadids) > $COMP(id) && $COMP(id) >= 0 && | if($PARENT(nbroadcastids) > $COMP(id) && $COMP(id) >= 0 | |||
i == $PARENT(threadids[$COMP(id)])) { | && | |||
i == $PARENT(broadcastids[$COMP(id)])) { | ||||
nthc += $COMP(whichdims_count); | nthc += $COMP(whichdims_count); | |||
} | } | |||
for(j=0; j<$COMP(whichdims_count); j++) { | for(j=0; j<$COMP(whichdims_count); j++) { | |||
if($COMP(whichdims[j] == i)) {flag=1; break;} | if($COMP(whichdims[j] == i)) {flag=1; break;} | |||
} | } | |||
if(flag) { | if(flag) { | |||
continue; | continue; | |||
} | } | |||
$CHILD(dims[nthc]) = $PARENT(dims[i]); | $CHILD(dims[nthc]) = $PARENT(dims[i]); | |||
$PRIV(incs[nthc]) = $PARENT(dimincs[i]); | $PRIV(incs[nthc]) = $PARENT(dimincs[i]); | |||
nthc++; | nthc++; | |||
} | } | |||
for(i=0; i<$COMP(whichdims_count); i++) { | for(i=0; i<$COMP(whichdims_count); i++) { | |||
int cdim,pdim; | int cdim,pdim; | |||
cdim = i + | cdim = i + | |||
($PARENT(nthreadids) > $COMP(id) && $COMP(id) >= 0? | ($PARENT(nbroadcastids) > $COMP(id) && $COMP(id) >= 0? | |||
$PARENT(threadids[$COMP(id)]) : $PARENT(ndims)) | $PARENT(broadcastids[$COMP(id)]) : $PARENT(ndims)) | |||
- $COMP(nrealwhichdims); | - $COMP(nrealwhichdims); | |||
pdim = $COMP(whichdims[i]); | pdim = $COMP(whichdims[i]); | |||
if(pdim == -1) { | if(pdim == -1) { | |||
$CHILD(dims[cdim]) = 1; | $CHILD(dims[cdim]) = 1; | |||
$PRIV(incs[cdim]) = 0; | $PRIV(incs[cdim]) = 0; | |||
} else { | } else { | |||
$CHILD(dims[cdim]) = $PARENT(dims[pdim]); | $CHILD(dims[cdim]) = $PARENT(dims[pdim]); | |||
$PRIV(incs[cdim]) = $PARENT(dimincs[pdim]); | $PRIV(incs[cdim]) = $PARENT(dimincs[pdim]); | |||
} | } | |||
} | } | |||
$SETDIMS(); | $SETDIMS(); | |||
PDL_RETERROR(PDL_err, PDL->reallocthreadids($CHILD_PTR(), | PDL_RETERROR(PDL_err, PDL->reallocbroadcastids($CHILD_PTR(), | |||
PDLMAX($COMP(id)+1, $PARENT(nthreadids)))); | PDLMAX($COMP(id)+1, $PARENT(nbroadcastids)))); | |||
for(i=0; i<$CHILD(nthreadids)-1; i++) { | for(i=0; i<$CHILD(nbroadcastids)-1; i++) { | |||
$CHILD(threadids[i]) = | $CHILD(broadcastids[i]) = | |||
($PARENT(nthreadids) > i ? | ($PARENT(nbroadcastids) > i ? | |||
$PARENT(threadids[i]) : $PARENT(ndims)) + | $PARENT(broadcastids[i]) : $PARENT(ndims)) + | |||
(i <= $COMP(id) ? - $COMP(nrealwhichdims) : | (i <= $COMP(id) ? - $COMP(nrealwhichdims) : | |||
$COMP(whichdims_count) - $COMP(nrealwhichdims)); | $COMP(whichdims_count) - $COMP(nrealwhichdims)); | |||
} | } | |||
$CHILD(threadids[$CHILD(nthreadids)-1]) = $CHILD(ndims); | $CHILD(broadcastids[$CHILD(nbroadcastids)-1]) = $CHILD(ndims); | |||
', | ', | |||
); | ); | |||
=head2 unthread | =head2 unbroadcast | |||
=cut | =cut | |||
pp_def( | pp_def( | |||
'unthread', | 'unbroadcast', | |||
Doc => <<'EOD', | Doc => <<'EOD', | |||
=for ref | =for ref | |||
All threaded dimensions are made real again. | All broadcasted dimensions are made real again. | |||
See [TBD Doc] for details and examples. | See [TBD Doc] for details and examples. | |||
=cut | =cut | |||
EOD | EOD | |||
P2Child => 1, | P2Child => 1, | |||
TwoWay => 1, | TwoWay => 1, | |||
AffinePriv => 1, | AffinePriv => 1, | |||
OtherPars => 'int atind;', | OtherPars => 'int atind;', | |||
RedoDims => ' | RedoDims => ' | |||
int i; | int i; | |||
$SETNDIMS($PARENT(ndims)); | $SETNDIMS($PARENT(ndims)); | |||
$DOPRIVALLOC(); | $DOPRIVALLOC(); | |||
$PRIV(offs) = 0; | $PRIV(offs) = 0; | |||
for(i=0; i<$PARENT(ndims); i++) { | for(i=0; i<$PARENT(ndims); i++) { | |||
int corc; | int corc; | |||
if(i<$COMP(atind)) { | if(i<$COMP(atind)) { | |||
corc = i; | corc = i; | |||
} else if(i < $PARENT(threadids[0])) { | } else if(i < $PARENT(broadcastids[0])) { | |||
corc = i + $PARENT(ndims)-$PARENT(threadids[0]); | corc = i + $PARENT(ndims)-$PARENT(broadcastids[0 | |||
]); | ||||
} else { | } else { | |||
corc = i - $PARENT(threadids[0]) + $COMP(atind); | corc = i - $PARENT(broadcastids[0]) + $COMP(atin d); | |||
} | } | |||
$CHILD(dims[corc]) = $PARENT(dims[i]); | $CHILD(dims[corc]) = $PARENT(dims[i]); | |||
$PRIV(incs[corc]) = $PARENT(dimincs[i]); | $PRIV(incs[corc]) = $PARENT(dimincs[i]); | |||
} | } | |||
$SETDIMS(); | $SETDIMS(); | |||
', | ', | |||
); | ); | |||
pp_add_exported('', 'dice dice_axis'); | pp_add_exported('', 'dice dice_axis'); | |||
pp_addpm(<<'EOD'); | pp_addpm(<<'EOD'); | |||
skipping to change at line 2330 | skipping to change at line 2338 | |||
} /* end of arg-parsing loop */ | } /* end of arg-parsing loop */ | |||
$COMP(idim_top) = idim; | $COMP(idim_top) = idim; | |||
$COMP(odim_top) = odim; | $COMP(odim_top) = odim; | |||
SLICE-MC | SLICE-MC | |||
, | , | |||
RedoDims => q{ | RedoDims => q{ | |||
PDL_Indx i; | PDL_Indx i; | |||
PDL_Indx PDIMS; | PDL_Indx PDIMS; | |||
int o_ndims_extra = PDLMAX(0, $PARENT(ndims) - $COMP(idim_top) ); | int o_ndims_extra = PDLMAX(0, $PARENT(ndims) - $COMP(idim_top) ); | |||
/* slurped dims from the arg parsing, plus any extra thread di ms */ | /* slurped dims from the arg parsing, plus any extra broadcast dims */ | |||
$SETNDIMS( $COMP(odim_top) + o_ndims_extra ); | $SETNDIMS( $COMP(odim_top) + o_ndims_extra ); | |||
$DOPRIVALLOC(); | $DOPRIVALLOC(); | |||
$PRIV(offs) = 0; /* Offset vector to start of slice */ | $PRIV(offs) = 0; /* Offset vector to start of slice */ | |||
for(i=0; i<$COMP(nargs); i++) { | for(i=0; i<$COMP(nargs); i++) { | |||
/** Belt-and-suspenders **/ | /** Belt-and-suspenders **/ | |||
if( ($COMP(idim[i]) < 0) && ($COMP(odim[i]) < 0) ) { | if( ($COMP(idim[i]) < 0) && ($COMP(odim[i]) < 0) ) { | |||
PDL->changed($CHILD_PTR(), PDL_PARENTDIMSCHANGED, 0); | PDL->changed($CHILD_PTR(), PDL_PARENTDIMSCHANGED, 0); | |||
$CROAK("Hmmm, both dummy and squished -- this can never happen. I quit."); | $CROAK("Hmmm, both dummy and squished -- this can never happen. I quit."); | |||
} | } | |||
skipping to change at line 2399 | skipping to change at line 2407 | |||
if(!inc) | if(!inc) | |||
inc = (start <= end) ? 1 : -1; | inc = (start <= end) ? 1 : -1; | |||
$CHILD( dims )[ $COMP(odim)[i] ] = PDLMAX(0, (end - s tart + inc) / inc); | $CHILD( dims )[ $COMP(odim)[i] ] = PDLMAX(0, (end - s tart + inc) / inc); | |||
$PRIV( incs )[ $COMP(odim)[i] ] = $PARENT(dimincs)[ $COMP(idim)[i] ] * inc; | $PRIV( incs )[ $COMP(odim)[i] ] = $PARENT(dimincs)[ $COMP(idim)[i] ] * inc; | |||
$PRIV(offs) += start * $PARENT( dimincs )[ $COMP(idim )[i] ]; | $PRIV(offs) += start * $PARENT( dimincs )[ $COMP(idim )[i] ]; | |||
} /* end of normal slice case */ | } /* end of normal slice case */ | |||
} /* end of trapped strange slice case */ | } /* end of trapped strange slice case */ | |||
} /* end of non-dummy slice case */ | } /* end of non-dummy slice case */ | |||
} /* end of nargs loop */ | } /* end of nargs loop */ | |||
/* Now fill in thread dimensions as needed. idim and odim per sist from the parsing loop */ | /* Now fill in broadcast dimensions as needed. idim and odim persist from the parsing loop */ | |||
/* up above. */ | /* up above. */ | |||
for(i=0; i<o_ndims_extra; i++) { | for(i=0; i<o_ndims_extra; i++) { | |||
$CHILD(dims)[ $COMP(odim_top) + i ] = $PARENT(dims)[ $COMP(i dim_top) + i ]; | $CHILD(dims)[ $COMP(odim_top) + i ] = $PARENT(dims)[ $COMP(i dim_top) + i ]; | |||
$PRIV(incs)[ $COMP(odim_top) + i ] = $PARENT(dimincs)[ $COMP (idim_top) + i ]; | $PRIV(incs)[ $COMP(odim_top) + i ] = $PARENT(dimincs)[ $COMP (idim_top) + i ]; | |||
} | } | |||
$SETDIMS(); | $SETDIMS(); | |||
} # end of RedoDims for slice | } # end of RedoDims for slice | |||
); | ); | |||
End of changes. 56 change blocks. | ||||
65 lines changed or deleted | 74 lines changed or added |