; file: crossalignaiacutoutfits2cubes.pro ; init: Nov 25 2013 Rob Rutten Sac Peak ; last: Nov 25 2013 Rob Rutten Sac Peak ; note: NO GOOD because the (x,y) location keywords are no good filepath='/media/usb0/alldata/SDO/2013-09-24-peter/level2/' selstringA='_1700' selstringB='_1600' outfileA='/media/usb0/alldata/SDO/2013-09-24-peter/cubes/newA.fits' outfileB='/media/usb0/alldata/SDO/2013-09-24-peter/cubes/newB.fits' ; channel A filesA=findfile(filepath+'*'+selstringA+'*') read_sdo,filesA[0],indA,imA nxA=indA.naxis1 nyA=indA.naxis2 ntA=n_elements(filesA) scaleA=indA.cdelt1 rotA=indA.crota2 xllA=indA.crpix1 yllA=indA.crpix2 ; channel B filesB=findfile(filepath+'*'+selstringB+'*') read_sdo,filesB[0],indB,imB nxB=indB.naxis1 nyB=indB.naxis2 ntB=n_elements(filesB) scaleB=indB.cdelt1 rotB=indB.crota2 xllB=indB.crpix1 yllB=indB.crpix2 nx=min([nxA,nxB]) ny=min([nyA,nyB]) nt=min([ntA,ntB]) ; set file type bigendian=1 bitpix=indA.bitpix ; define output header mkhdr,outheader,abs(bitpix)/8,[nx,ny,nt] sizeoutheader=size(outheader) ; fits outheader = Nx36 "card images" = Nx2880 bytes outheadersize=(1+fix(sizeoutheader[1]/36.))*2880 ; open output file A for assoc, write header get_lun, unit_outA if (bigendian) then openw,unit_outA,outfileA,/swap_if_little_endian $ else openw,unit_outA,outfileA if (bitpix eq -32) then outassocA=assoc(unit_outA,fltarr(nx,ny),outheadersize) if (bitpix eq 16) then outassocA=assoc(unit_outA,intarr(nx,ny),outheadersize) if (bitpix eq 8) then outassocA=assoc(unit_outA,bytarr(nx,ny),outheadersize) if (outheadersize ne 0) then begin rec=assoc(unit_outA, bytarr(outheadersize)) rec[0]=byte(outheader) endif ; open output file B for assoc, write header get_lun, unit_outB if (bigendian) then openw,unit_outB,outfileB,/swap_if_little_endian $ else openw,unit_outB,outfileB if (bitpix eq -32) then outassocB=assoc(unit_outB,fltarr(nx,ny),outheadersize) if (bitpix eq 16) then outassocB=assoc(unit_outB,intarr(nx,ny),outheadersize) if (bitpix eq 8) then outassocB=assoc(unit_outB,bytarr(nx,ny),outheadersize) if (outheadersize ne 0) then begin rec=assoc(unit_outB, bytarr(outheadersize)) rec[0]=byte(outheader) endif shiftx=xllB-xllA ; 4.61 ; pos is to right ;; much too large! shifty=yllB-yllA ; 0.69 xscale=scaleB/scaleA ; 0.99424864 YES! yscale=xscale rotation=rotB-rotA : -0.0006620 ; trials one by one ;; shiftx=0 ;; shifty=0 ;; xscale=0 ;; yscale=0 ;; rotation=0 ; doesn't do anything ;; for it=0,nt-1 do begin for it=0,2 do begin read_sdo,filesA[it],indA,imA read_sdo,filesB[it],indB,imB newimA=imA[0:nx-1,0:ny-1] newimB=imB[0:nx-1,0:ny-1] newimB=shift_img(newimB,[shiftx,shifty],$ xscale=xscale,yscale=yscale,rot=rotation,anchor=[xllA,-yllA]) outassocA[it]=newimA outassocB[it]=newimB endfor ; free the output files free_lun,unit_outA,unit_outB print,' === wrote ',outfileA,' and ',outfileB end