;Arranges a stack such that the corresponding pixels from all layers are in the same column
;  input
;    bbox - bounding box of the entire orbit. From the bbox field of the level 1B structure
;    scene_bboxes - bounding boxes of each scene. From the scene_bboxes field of the level 1B structure
;    data - array of pointers to layers. From the level 1B structure (albedo, scattering_angle, etc)
;  returns
;    A 3D array where each layer (result[*,*,i]) is data from one layer, and each column (result[x,y,*])
;    is data covering the same chunk of sky.
function arrange_stack,bbox,scene_bboxes,data,l1b
  ;Get the total number of layers from the size of the scene_bboxes array
  s=size(scene_bboxes,/dimensions)
  n_scenes=s[0]
  ;Set up the arranged
  stack=fltarr(bbox[2]-bbox[0]+1,bbox[3]-bbox[1]+1,n_scenes)*!values.f_nan

  alb_p=dereference_single(l1b.albedo)

  ;For each layer
  for i=0,n_scenes-1 do begin
    ;Find out where this layer fits into the whole swath
    scene_x=scene_bboxes[i,0]-bbox[0]
    scene_y=scene_bboxes[i,1]-bbox[1]
    scene_w=scene_bboxes[i,2]-scene_bboxes[i,0]+1
    scene_h=scene_bboxes[i,3]-scene_bboxes[i,1]+1

    ;put the layer data in its place
    if ptr_valid(data[i]) then begin       
       stack[scene_x:scene_x+scene_w-1,scene_y:scene_y+scene_h-1,i]=*(data[i])
    endif else begin
       albedo_mask = *alb_p[i]
       albedo_mask = albedo_mask/albedo_mask  ; = 1 where we have an albedo and NaN where we don't
       
       if strmatch(data[i],'px') then begin
          stack[scene_x:scene_x+scene_w-1,scene_y:scene_y+scene_h-1,i]=0.0*albedo_mask
       endif else if strmatch(data[i],'mx') then begin
          stack[scene_x:scene_x+scene_w-1,scene_y:scene_y+scene_h-1,i]=1.0*albedo_mask
       endif else if strmatch(data[i],'py') then begin 
          stack[scene_x:scene_x+scene_w-1,scene_y:scene_y+scene_h-1,i]=2.0*albedo_mask
       endif else if strmatch(data[i],'my') then begin 
          stack[scene_x:scene_x+scene_w-1,scene_y:scene_y+scene_h-1,i]=3.0*albedo_mask
       endif else begin
          stack[scene_x:scene_x+scene_w-1,scene_y:scene_y+scene_h-1,i]=data[i]*albedo_mask
       endelse
    end

  end
  return,stack
end