Re: A little trouble with very large arrays.

[Niklas Holsti <niklas.holsti@tidorum.invalid>]

On 18-10-05 19:56 , Shark8 wrote:
> On Friday, October 5, 2018 at 12:36:47 AM UTC-6, Dmitry A. Kazakov wrote:
>> On 2018-10-04 23:38, Shark8 wrote:
>>> GNAT won't even compile an array like this [999 indexes].
>>>
>>> What's the proper way to go about doing this?
>>
>> A wrong way dealing with protocols is attempting to define an Ada type
>> having the exact representation of the data as defined by the protocol.
>> It is both useless and difficult to impossible, especially if bits are
>> involved.
>
> Protocol?
> FITS is a file-format.

Which is defined as a sequence of bytes, so a FITS file is equivalent to 
a message in a protocol such as SMTP etc. Usually a very _long_ message, 
of course.

> The only reason bits are involved at all in the
> spec is because it was developed back when some machines had 9-bit
> bytes.

FITS version 4.0 defines everything with 8-bit bytes, as far as I could 
see with a glance at the standard. Do you need to process some older 
FITS files with a different byte-size?

Yes, the FITS block size (2880 octets) was chosen to be divisible by 9, 
and other ancient word-sizes and byte-sizes, but so what?

>> You could use a flat array internally and provide operations for image
>> serialization/deserialization in whatever format, e.g. by
>> Get_Pixel/Set_Pixel.
>
> I tried this, it doesn't quite work though. (Stack overflow, oddly enough.)
>     Function Flatten( Item : Axis_Dimensions ) return Natural is
>       (case Item'Length is
> 	   when 0 => 1,
> 	   when 1 => Item( Item'First ),
> 	   when 2 => Item( Item'First ) * Item( Item'Last ),
> 	   when others =>
> 	     Flatten( Item(Item'First..Item'Last/2) ) *

That Item'Last/2 does not seem right. If you want the middle index, it 
should be (Item'First + Item'Last) / 2. Perhaps this error leads to an 
unending recursion, explaining the stack overflow.

> 	     Flatten( Item(Axis_Count'Succ(Item'Last/2)..Item'Last) )
>       );

But what is this function supposed to do? Is it meant to compute the 
length (number of elements) in the flattened array? That is just the 
product of the Axis_Dimension values, isn't it?

    function Product (Item : Axis_Dimensions) return Natural
    is
       Result : Natural := 1;
    begin
       for I in Item'Range loop
          Result := Result * Item(I);
       end loop;
       return Result;
    end Product;

For computing the position (flattened index) of an element in a 
flattened multi-dimensional array, you need a function that takes two 
arguments:
- a vector giving the length of each axis
- a vector giving the index (of the element) for each axis.

Coding that function as a double recursion gives no benefit IMO. A 
simple loop is better, as in the function above.

Also remember that the FITS array is in Fortran order, so the index of 
the first axis varies most rapidly in the flattened sequence of array 
elements. This can be done by a "loop .. in reverse ...".

-- 
Niklas Holsti
Tidorum Ltd
niklas holsti tidorum fi
       .      @       .