subr overlap_diff, im1, im2, xc1, yc1, xc2, yc2, im3, xc3, yc3, sc
 ;returns a difference image of the overlapping area, im3 = im1 - im2
 ;im3 will be dimensioned to the overlap area and xc3 and yc3
 ;are its center on the sun (or whatever was used for xc1, etc)
 ;only nearest neighbor alignment is done, no fractional pixel
 ;shifts, sc is the scale between pixels and the coord. system, usually
 ;xc1, etc are in arc seconds which means that sc will be "/pixel, both
 ;images must have the same scale
 ;if there is no overlap a 1x1 image is returned with a position that
 ;is the average of the two input - this is an attempt to make it
 ;sensible for a series that has a few non overlapping cases
 ;should be safe for im3 to be the same name as either im2 or im1
  
 ;dump, im1, im2
 nx1 = dimen(im1,0)	ny1 = dimen(im1,1)
 nx2 = dimen(im2,0)	ny2 = dimen(im2,1)
 if sc le 0.0 then sc = 1.0
 scinverse = 1./sc
 xq = (xc1 - xc2) * scinverse
 nqx = 0.5*(nx1 - nx2)
 yq = (yc1 - yc2) * scinverse
 ;11/22/2000 - bug in next line had ny1 and ny2 reversed
 nqy = 0.5*(ny1 - ny2)
 if abs(xq) gt (nqx + nx2 + 1) or abs(yq) gt (nqy + ny2 + 1) then {
  ty,'no overlap in overlap_diff'
  xc3 = .5*(xc1 + xc2)		yc3 = .5*(yc1 + yc2)
  im3 = zero(intarr(1,1))
  return
 }
 xq = xq - nqx
 if xq ge 0 then {
 	x0_1 = 0   x0_2 = rfix(xq)  x1_2 = nx2 - 1  x1_1 = x1_2 - x0_2
 	if x1_1 ge nx1 then { x1_1 = nx1 - 1   x1_2 = x1_1 + x0_2 }
 } else {
 	x0_2 = 0   x0_1 = rfix(-xq) x1_1 = nx1 - 1  x1_2 = x1_1 - x0_1
 	;but if 2 contained within 1, then that could go out of range
 	if x1_2 ge nx2 then { x1_2 = nx2 - 1   x1_1 = x1_2 + x0_1 }
 }
 yq = yq - nqy
 if yq ge 0 then {
 	y0_1 = 0   y0_2 = rfix(yq)  y1_2 = ny2 - 1  y1_1 = y1_2 - y0_2
 	if y1_1 ge ny1 then { y1_1 = ny1 - 1   y1_2 = y1_1 + y0_2 }
 } else {
 	y0_2 = 0   y0_1 = rfix(-yq) y1_1 = ny1 - 1  y1_2 = y1_1 - y0_1
 	;but if 2 contained within 1, then that could go out of range
 	if y1_2 ge ny2 then { y1_2 = ny2 - 1   y1_1 = y1_2 + y0_1 }
 }
 
 ;compute the center of the subarea in the global coordinates
 xc3 = xc1 + sc*0.5*(x0_1 + x1_1 + 1 - nx1)
 yc3 = yc1 + sc*0.5*(y0_1 + y1_1 + 1 - ny1)
 
 ;extract the subimages, always use at least I*2 here, only need to apply
 ;to one array because of auto upgrade
 ;;ty,'new center = ', xc3, yc3
 ;;ty,'cutouts for 1:',x0_1,x1_1, y0_1,y1_1
 ;;ty,'cutouts for 2:',x0_2,x1_2, y0_2,y1_2
 xq = im2(x0_2:x1_2, y0_2:y1_2)
 ;;tv,xq,0,0,$movie_show_win
 if symdtype(xq) eq 0 then xq = word(xq)
 ;;tv,im1(x0_1:x1_1, y0_1:y1_1),0,0,31
 xq = im1(x0_1:x1_1, y0_1:y1_1) - xq
 switch, im3, xq
 ;;ty,'output array'
 ;;dump, im3
 xq = 0
 return
 endsubr
 ;========================================================