var('x,y')
def region(x,y):
if x^2 + y^2 < 16*y and x^2 + y^2 > 64:
return 1
else:
return 0
(contour_plot(region,(x,-9,9),(y,-9,18),cmap=[(1,1,1),(1,0,0)],plot_points=300) +\
implicit_plot(x^2+y^2==16*y,(x,-9,9),(y,-9,18),cmap=[(0,0,1)],linewidth=2) +\
implicit_plot(x^2+y^2==64,(x,-9,9),(y,-9,18),linewidth=2))\
.show(aspect_ratio=1,frame=False,axes=True)
